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Charm final states in deep inelastic scattering constitute ~ 25% of the inclusive cross-section at small 
x as measured at HERA. These data can reveal important information on the charm and gluon structure 
of the nucleon if they are interpreted in a consistent perturbative QCD framework which is valid over the 
entire energy range from threshold to the high energy limit. We describe in detail how this can be carried 
out order-by-order in PQCD in the generalized MS formalism of Collins (generally known as the ACOT 
approach), and demonstrate the inherent smooth transition from the 3-flavor to the 4-flavor scheme in a 
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1 Introduction 



Recent measurements of charm production in deep inelastic scattering (DIS) at HERA |l], ^ have shown 
that up to 25% of the total cross-section at small-x contains charm in the final state. This is within the 
expectations of perturbative QCD based on conventional parton distributions. We are now in a position to 
utilize this process to study details of the production mechanism of heavy quarks in general, and to extract 
useful information on the charm and gluon structure of the proton in particular. 

Conventional perturbative QCD (PQCD) theory is formulated in terms of zero-mass quark-partons. For 
processes depending on one hard scale Q, the well-known factorization theorem then provides a straightfor- 
ward procedure for order-by-order perturbative calculations, as well as an associated intuitive parton picture 
interpretation of the perturbation series. Heavy quark production presents a challenge in PQCD because 
the heavy quark mass, mH (H = c,b,t), provides an additional hard scale which complicates the situation 
- it requires a different organization of the perturbative series depending on the relative magnitudes of mH 
and Q. 

The two standard methods for PQCD calculation of heavy quark processes represent two diametrically 
opposite ways of reducing the two-scale problem to a one-scale problem, (i) In the conventional parton model 
approach used in many global QCD analyses of parton distributions j| [I| ||] and Monte Carlo programs, the 
zero-mass parton approximation is applied to a heavy quark calculation as soon as the typical energy scale^] 
of the physical process Q is above the mass threshold mH. This leaves Q as the only apparent hard scale 
in the problem, (ii) In the heavy quark approach which played a dominant role in "NLO calculations" of 
the production of heavy quarks || [7], || , the quark H is always treated as a "heavy" particle and never as a 
parton. The mass parameter mH is explicitly kept along with Q as if they are of the same order, irrespective 
of their real relative magnitudes. 

The co-existence of these two opposite approaches represents an uneasy dichotomy in the current liter- 
ature. On physical grounds, the zero-mass parton picture of heavy quarks should be applicable at energy 
scales very much larger than the relevant quark mass, mH <C Q, whereas the heavy quark approach (often 
referred to as the fixed-flavor-number (FFN) scheme) should be more appropriate at energy scales compa- 
rable to the quark mass mH ~ Q. The actual experimental regime often lies in between these two extreme 
regions, where the validity of cither approach can be called into question. There is, however, a natural way to 
incorporate both approaches in a unified framework in PQCD which provides a smooth transition between 
the two. This has been formulated in a series of papers over the years, Refs. ||, [H], [if], [l^, which 
has now been adopted, in different guises, by most recent literature on heavy quark production in PQCD. 



To see the basic ideas behind this unified picture, let us focus explicitly on the production of charm 
(H = c) in deep inelastic scattering. All considerations apply to a generic heavy quark. Consider the PQCD 
calculation of the F2(x, Q) structure function which receives substantial contribution from charm production 
as mentioned earlier. The underlying physical ideas are illustrated graphically in Fig. ^ where the charm 
contribution to this structure function, denoted F2 c (x,Q), is plotted as a function of Q at some fixed value 
of x. Near threshold Q ~ mc, it is natural to consider the charm quark as a heavy particle, and to adopt the 
3 active-parton- flavor scheme of calculation (the "heavy quark approach" ) . As Q becomes large compared to 
mc, the fixed-order calculation in this approach becomes unreliable since the perturbative expansion contains 
terms of the form as n log n (mc 2 /Q 2 ) at any order n, which ruin the convergence of the series — these terms 
are not infra-red safe as mc — > or Q — > oo. Thus the uncertainty of the 3-flavor calculation grows as 
Q/mc becomes large. This is illustrated in Fig. [l] as an error band marked by horizontal hashes which is 
narrow near threshold but becomes ever wider as Q/mc increases. On the other hand, starting from the high 
energy end (Q 3> mc), the most natural calculational scheme to adopt is the conventional 4-flavor scheme 
with active charm partons. (In this approach, the infra-red unsafe large logarithms mentioned earlier are 
"resummed" and absorbed into the finite charm parton distributions.) However, as we go down the energy 
scale toward the charm production threshold region, the 4-flavor calculation becomes unreliable because the 
approximation mc = deteriorates as Q — > mc. The uncertainty band of such a calculation is outlined in 
Fig. by the vertical hashes - it is narrow at high energies, but becomes increasingly wider as one approaches 
the threshold region. 

d We use Q as the generic name for a typical kinematic physical scale. It could be Q, W, or pT, depending on the process. 
6 See Ref. Jw|] for a brief review. 
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Figure 1: Schematic illustration of the uncertainty bands (hence the relative merits) of the 3-fiavor and 
4-flavor scheme calculations of the charm contribution to the inclusive structure function - and how they 
can be naturally merged into a composite scheme. 



The intuitive ideas embodied in Fig [l] illustrate that: (i) these two conventional approaches are individ- 
ually unsatisfactory over the full energy range, but are mutually complementary; and (ii) the most reliable 
PQCD prediction for the physical F2(x,Q), at a given order of calculation, can be obtained by utilizing 
the most appropriate scheme at that energy scale Q, resulting in a composite scheme, as represented by 
the cross-hashed region in Fig. § The use of a composite scheme consisting of different numbers of flavors 
in different energy ranges, rather than a fixed number of flavors, is familiar in the conventional zero-mass 
parton picture. The new formalism espoused in Refs. |], |h], [Uj [l2| provides a quantum field theoretical 
basis JlQj , p"3[ for this intuitive picture in the presence of non-zero quark mass. 

This generalization brings about several important distinguishing features and insights. First, after 
the hard cross-section is rendered infra-red safe by factorizing out the ln(mc) terms, the remaining finite 
dependence on mc can be kept in the hard cross-sections to maintain better accuracy in the intermediate 
energy region.^One can then show that the unified formalism reproduces the two conventional approaches in 
well-defined ways in their respective regions of applicability @-(l3). Secondly, it has become clear recently 
that there is much inherent flexibility in the choice of the transition energy scale (cf. Fig. 0) as well as 
the detailed matching condition between the 3- and 4-flavor calculations. This makes it possible to have a 
variety of different implementations of the new formalism fl5| , [u| 18 , with different emphases. Properly 
understood, this feature can be exploited to lend strength to the formalism; but, during the development of 
the theory, the subtleties have given rise to much misunderstanding and confusion in existing literature. 

In this paper, we study charm production in DIS using the general mass formalism, and compare the 
results with recent HERA data. The main goals of this work are: 

(i) Through this specific example, we give a concise and careful presentation of the general formalism, 
including the relatively simple operational procedure for calculating its various components at any order of 
a s . We hope this will help fill the gap between the original, relatively sketchy, ACOT paper || and the 
recent, more technical, all-order proof of the formalism by Collins Jig] . 

(ii) We carry out the numerical calculations to make concrete the intuitive ideas illustrated in Fig. Q and to 
demonstrate the validity of the physical principles underlying the composite scheme. 

(iii) We show that the flexibility of the general formalism mentioned above can result in efficient PQCD 
calculations of inclusive quantities at relatively low order in a s compared to FFN calculations - because the 
relevant physics has been effectively "resummed" by the appropriate scheme adopted for the given energy 
scale. 

Item (i) is basic; it forms the foundation for the other two parts. However, since this part is about 
the clarification of the existing theoretical formalism rather than new work, and since not all readers are 
concerned with theoretical precision, we have elected to place it in the Appendix. Hopefully, this will increase 

f By contrast, in the standard literature the mc — > limit is tied to the proof of factorization in the first place. This 
association is not needed, cf. Jli| 
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the accessibility of the main body of the paper. In regards to the underlying physics, we believe the discussion 
in the appendix should make a useful contribution to clear the confusion and misunderstanding among the 
various approaches that have been proposed in the recent literature following @ . 

Based on the terminologies discussed in this Appendix, in Sec. || we describe the complete order a s 
calculation carried out in this paper in relation to previous work on this subject. The new calculation 
extends the validity of the original ACOT results to NLO in the high energy regime - on the same level as the 
conventional zero-mass total inclusive structure functions. In addition, the new perspective, as illustrated 
in Fig. [I] above and discussed quantitatively in the paper proper, allows a re-assessment of the physical 
predictions of the order a s calculation near the threshold region, making it a viable alternative to the 
order a 2 FFN calculation. In Sec. |^ we present the numerical results on inclusive charm production, and 
demonstrate that the validity and the efficiency of the general formalism, as described in (ii) and (iii) above, 
are indeed seen at this order. We show that very good agreement with recent HERA data on inclusive 
F2 C is obtained in practice. For this study, we have developed a new implementation of the generalized 
MS formalism using Monte Carlo methods. This implementation allows the computation of differential 
distributions, with kinematic cuts, such as da D /dpt, da D /dQ, and da D jdr\ which we present in Sec. || 
These results, are in qualitative agreement with available data from HERA. However, they are not as good as 
those of the order as 2 calculation in the 3-fl.avor scheme. This is to be expected for differential distributions 
with experimental cuts, since the (resummed) low-order calculation contains more severe approximations to 
the kinematics of the final state partons. For future quantitative studies, the general formalism needs to 
be expanded to incorporate higher-order results (adaptable from existing FFN calculations). This point is 
discussed in the concluding section, along with other observations. 



2 Total Inclusive Structure Functions in the general formalism 

We consider the inclusive DIS structure functions, such as F 2 , focusing on the contribution of a massive 
quark. For definiteness, we assume that the only relevant quark with non-zero mass is the charm quark. 
The generic leptoproduction process is depicted in Fig. ||: 

h + A^e 2 + x, (i) 

where A is a hadron, are leptons, and X represents the summed-over final state hadronic particles. Note 
that X may or may not contain a visible heavy-flavor hadron. After the calculable leptonic part of the 




Figure 2: The generic inclusive lepto-production process. We are primarily interesting in final states which 
contain some heavy quarks, in particular charm. 

cross-section has been factored out, we work with the hadronic process induced by the virtual vector boson 
7* of momentum q and polarization A: 

1 *(q,X)+A(P) -^X(P X ). (2) 
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Although our considerations apply to DIS processes induced by W and Z as well, we shall explicitly refer 
to the neutral current interaction with the exchange of a virtual photon 7* in order to be concrete. The 
cross-section is expressed in terms of the hadronic tensor 

Wxa(q,P) = ^ £ (P\e x -J ] \Px){2nf5 {i) {P + q-Px)(Px\e a -J\P). (3) 

X(Px),spin 

where denotes a sum over all final hadronic states. In most cases, it suffices to consider the diagonal 
elements of the tensor F x = W\\(q, P). 

The factorization theorem in the presence of non-zero quark masses - assumed in jl2| and established to 
all orders in PQCD - states that the inclusive cross-section can be written as a convolution: 

Fi(Q 2 , x, ..) =J2 f%(x, fi) ® £ q ,a(x, Q/fM, Q/m c , a s (fx)) + o(k/Q) P (4) 

a 

where /j| is the distribution of parton a inside the hadron A, uj a .x is the perturbatively calculable hard 
cross-section for 7* + a — » X, p is some positive number, fi denotes collectively the renormalization and 
factorization scales, and a convolution in the x variable is implied. The helicity structure functions F x are 



simply related to the familiar -Fi,2,3 H 



The exact way that the physical structure function factorizes into the long-distance (/^) and the short- 
distance (fi a ,x) pieces on the right-hand side of Eq. (Q) depends on the scheme used to define the parton 
distributions. The physical structure function should be independent of any calculational scheme; there- 
fore, the definition of the hard cross-sections is determined by the subtraction procedure used to define the 
parton distributions. As discussed in the Introduction, the general formalism consists of the 3-flavor scheme 
at low energy scales, the 4-flavor scheme at high energy scales, and a suitably chosen transition region where 
matching conditions between the two schemes are applied. The precise description of these elements of the 
formalism is given in the Appendix. Operational definitions of quantities needed in subsequent discussions 
are also discussed in some detail there. 



2.1 Previous calculations 

To put the current calculation in context, we first summarize the existing calculations of leptoproduction of 
charm using the precise definitions given in the Appendix. 

• NLO 3-flavor (3a 2 ) calculation [fr|: Most dedicated calculations of heavy flavor production in recent 
years have been carried out in this scheme. The LO process is 0(a\) heavy- flavor creation (HC), j*g — > cc. 
The NLO processes consist of the 0(a 2 s ) virtual corrections to 7*17 — ► cc as well as the real HC 7*/ — > ccl 
process, where I denotes any light parton. Cf. Fig. ||. This calculation becomes questionable when Q 3> m c 
- i.e., it ceases to be "NLO" in accuracy, as indicated in Fig. 0, because the perturbative expansion is 
actually in a s ln(Q 2 /ro 2 .) for large Q. 

• Zero-mass 4-flavor (ZM) (4a*; m c = 0) calculation: This is the formalism used in most conventional 
QCD parton model calculations and popular Monte Carlo programs. It represents an approximation to the 

general mass (GM) 4-flavor scheme by setting m c = in the hard cross-sections a) Q (x, Q, ^jN/-*) m '^> 
w™ c=0 (a;, ^jM)- The LO contribution consists of the O(a s ) j*c — > c heavy-flavor excitation (HE) process. 
The NLO contribution consists of the 0(a\) virtual corrections to j*c — > c plus the real HC j*g — > cc and, 
7*c — > gc HE processes. Cf. Fig. ^. This calculation is unreliable in the threshold region, as mentioned in 
the introduction. 

• LO generalized MS calculation (ACOT) [^2|: This represents the simplest implementation of the 
generalized formalism It emphasizes the overlapping physics underlying the 3-flavor and 4-flavor 
calculations in the region not far above threshold. The 4-flavor calculation consists of HE 7*c — > c 
(Fig. [|a) plus m c ^ HC 7*17 — > cc (Fig. ||d), with the requisite subtraction term which removes the 
m c mass-logarithm. Mathematically, the result on F 2 (x,Q, fi) can be shown to match that of the 0(a\) 
3-flavor scheme calculation (HC 7*3 — > cc) as the (unphysical) factorization scale fi approaches m c from 
above [[l2| . Physically, the predicted behavior of F^ for Q ~ m c will depend on the choice of /1 as a 
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Figure 3: Representative Feynman diagrams contributing to the calculation of the partonic structure func- 
tions for charm production in the 3-fLavor scheme due to: (a) order a\ heavy-flavor creation (HC) or gluon 
fusion mechanism; (b) order a s virtual correction to (a); and (c) order a 2 s real corrections involving gluon 
and light-quark initial states. 





Figure 4: Representative Feynman diagrams contributing to the calculation of the partonic structure func- 
tions for charm production in the 4-fLavor scheme due to: (a) order a° s heavy-flavor excitation (HE) mech- 
anism; (b) order a], virtual correction to (a); (c) order a], real corrections to the HE mechanism; and (d) 
order a), heavy-flavor creation mechanism. 



function of the physical variables, cf. next section. At a high energy scale, [i ~ Q 3> m c , this calculation 
only approximates the 4-fiavor NLO results: it contains the most important 0{a].) term, HC j*g — » cc 
(because of the large gluon distribution and the need for matching), but does not include the smaller 0(a*) 
terms represented by Fig. [|b,c. This calculation has been further studied by Kretzer and Schienbein pl| 
and Kramer et al [ |22| . 

• Variations on the "variable flavor number" theme: In recent years, the general approach proposed 
in [[)[ [l^] has been adopted by other groups, starting from different historical perspectives, Jl5|, |l6[ [Tt], 
In contrast to the fixed flavor-number (FFN) approach, these are usually referred to as being in the variable 
flavor-number (VFN) scheme. This terminology has caused some confusion, since although the common 
theme is that of ||, ^2|, as shown in Fig. [j], the various implementations differ considerably. Whether 
a particular implementation is self-consistent within the general framework of PQCD, or whether two 
implementations are compatible within the accuracy of the perturbative expansion, are often unclear or 
controversial (e.g. (23)) because of the complexity of the multi-scale problem and because of possible 
misunderstandings. It is beyond the scope of this paper to critically review these approaches. By pursuing 
the goals described in the introduction, we hope to provide a clearer picture of the general formalism of 
fl5}| , hence a better basis to help address some of the controversial issues in the future. 

2.2 The full order a s generalized MS calculation 

The calculation of charm contribution to the inclusive structure functions reported in this paper completes 
the order a s calculation in the general formalism, initiated in jl2| , including all the hard processes of Fig. [|. 
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The additional terms, although relatively small numerically at current energies, are required to make the 
4-flavor part of this calculation truly NLO at high energies, so that it becomes equivalent to the conventional 
zero-mass NLO theory used in most modern analyses of precision DIS data. In the remainder of this section 
we discuss the theoretical issues and uncertainties in the order a s 4-flavor calculation at all energy scales, in 
order to address the issues highlighted at the end of the introduction, Sec. [|. 

At order a\ the 3-flavor component of the composite scheme consists of only the (unsubtracted) 0(ct\) 
HC 7*5 — > cc process. The result is standard. Therefore, our main calculation concerns the J^-flavor scheme 
component. At order a\ the right-hand side of Eq. (Q) consists of three terms 

F^(Q 2 ,x,m c ,/i) = /1®°£ C C A 

+ \% (5) 

+ light — parton contributions , 

where the superscript (4) indicates that this is a 4-flavor calculation, and the hard- scattering cross- sections 
f° r the various subprocesses are calculated from the corresponding partonic cross- sections 1 <*)q\ (with- 
out the hat) according to the procedures described in the Appendix. A description of each of the terms 
follows: 

• The leading order 7*c — > c (Fig. ^a) partonic cross-section uj^, A is infra-red safe, thus 



c,X — w c,X ■ 



(6) 



• The 7*g — > cc (Fig. ^d) partonic cross-section contains a single power of In /i 2 /m 2 which can be factorized 
into the charm distribution function by the subtraction jl2j 



i s% x = i < A - v;® x,a, (7) 

where the cancelling logarithm with mass-singularity resides in the 0(a s ) perturbative parton distribution 
function 

l f c g = {a s /2^)P g ^ q {x)H^/rnl) (8) 

The virtual correction to 7*c — > c (Fig. ||b) plus the real 7*c — > gc (Fig. ||c) partonic process also contain 
In (/i 2 /m 2 ) terms which are factorized into the charm distribution function by the subtraction 

1 Z c c X x= ^l- 1 /^ ^, (9) 



where the logarithm appears in the 0{a 1 s ) perturbative parton distribution function^] 

(10) 



lrc = ^£.4 

Jc 2vr3 



1 +x2 ^(\ n ^4-l-2In(l-ar) 



1 — x I \ m 



Substituting Eqs. (|6|)-(^|) into Eq. (||), the right-hand side can be re-organized as 



+ (fi a® v s c -f%® v c c ) ® x,a (in 

, rr. ^ 1 rX \ > 



-r/1® X "l X X 
+ light — quark terms , 

where the In (/i/m c ) terms in the Lo a ,x factors are kept intact, and the needed subtraction terms are explicitly 
grouped with the leading 2— >1 term with the same kinematics. This is the form we use for the actual numerical 
calculations, which we implement using a Monte Carlo approach. 

s We have calculated these terms keeping a finite charm quark mass. As discussed in Ref. jljj, it would also have been 
consistent to calculate diagrams with an initial state charm quark using a zero charm quark mass. The errors near threshold 
in both methods of calculation are comparable and of order . 
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It is useful to compare this calculation with the NLO 3-flavor calculation Q which can be written, in the 
same notation, as follows 

Ff x {Q\x,m c ^) = fl® 1 ^ 

+ n ® 2 <l (i 2 ) 

+ light — quark terms 

The term common to the two schemes is O (a£) 7*<7 — > cc, which appears as the first line in both Eq. [ll] and 
Eq. [l^. For this comparison, one can consider the rest of the terms in these equations as complementary 
"corrections" to the first term. In particular, the last two terms in the 3-flavor formula Eq. [l^ are genuine 
O (oq) corrections to the common term in the threshold region; hence are commonly referred to as NLO. 
But at high Q 2 3> m 2 , these terms contain large logarithm factors In Q 2 /m 2 which vitiates the perturbation 
expansion. They are not NLO in this region. 

In the 4- flavor calculation as organized in the form Eq. [ll], the last two lines are also effectively O (a 2 ). 
This is because the distribution j c A is effectively O (a*) near threshold, and there is a built-in cancellation 
between the leading order 7*c — > c partonic cross-section and the first subtraction term in Eq. [ll], as discussed 
in detail in Ref.jjJ]]. Although these O (of) "correction" terms to the O (a*) common term (first line) do 
not contain the full NLO O (a 2 ) corrections at threshold, they do contain all the O (o^) contributions which 
are enhanced by In (Q 2 /m 2 ) and quickly dominate as Q 2 increases. In fact, these large logarithmic terms 
have been resummed to all orders in perturbation theory via the DGLAP-evolved charm parton distribution 
contributions in the second and third lines. Therefore, in the large Q 2 region the 4-flavor result represents 
a true NLO calculation. 

In principle, if carried out to all orders in a s , the 3-flavor and 4-flavor calculations in the form of 
F^\(Q 2 , x, m c , pi), i = 3, 4, would give exactly the same prediction for all values of the arguments. The de- 
pendence on the scheme (i) and the scale (p) arises from the truncation of the perturbation series. Therefore, 
in order to produce a physical prediction F^ h ^ s (Q 2 , x, m c ) and to relate the predictions in the two schemes, 
a number of additional steps must be taken. In the presence of a non-zero heavy quark mass, some of these 
steps are non-obvious; hence they can be the source of confusion. An explicit discussion of these elements 
of the calculation will make clear the flexibility as well as the uncertainties inherent in the formalism. This 
we do in the Appendix, as part of the more precise description of the formalism. Here, we mention only two 
features which are particularly relevant for the subsequent discussions of physical predictions. 

First, within each scheme (i = 3 or 4), one needs to specify \i as a function of the physical variables in 
order to make a physical prediction, i.e. 

F^ s {Q 2 ,x,m c ) = F ( ^ X (Q 2 ,x,m c , fj,(x,Q,m c )) 

Although there is considerable freedom in choosing fi(x, Q, m c ), two conditions should be met so that the pre- 
diction can be reliable: (i) /i must be of the order of Q or m c so that PQCD applies, and (ii) F A \(Q 2 ,x, m c , fi) 
must be relatively stable with respect to variations of /x for the (x, Q)-range of interest. This is the well-known 
scale- dependence of any PQCD calculation. For the problem at hand, a common choice for fi(x, Q, m c ) is 
\/ Q 2 + m 2 : it represents the typical virtuality of the internal parton lines in the important subprocesses. 
The presence of the uncertainty associated with the choice of ii(x, Q, m c ) in each scheme is illustrated in 
Fig. [I] by the respective bands .Fl 

Secondly, as shown in Fig. M one needs to identify an appropriate scale at which the 3-flavor and the 
4-flavor scheme predictions are both reasonable and mutually comparable, so that the transition from one to 
the other in the composite scheme can be made smoothly. In the next section we show that, for the inclusive 
charm production cross-section, these conditions can be met over a rather large range of Q, extending down 

h In the original ACOT paper, a different scale function fi(x,Q,m c ) was chosen: it was designed to enforce the condition 
_F( 4 ) (Q 2 , x, m c , fi(x, Q, m c )) — > F^q (Q 2 , x,m c ,m c ) as Q — > m c . In retrospect, this choice is artificial and unnecessary. The 
general formalism automatically ensures that, F^{Q 2 ,x,m c ,fi) -> {Q 2 , x , m c , m c ) as fi -> m c for given (Q 2 ,x) to the 

order which we are working; but it does not place any restriction on the behavior of F phys (Q 2 , x,m c ) as Q — > m c in any 

given scheme. As shown in Sec. ^, the more natural choice of scale fj,(x,Q,m c ) = */ Q 2 + "i- 2 - leads to much improved physical 
predictions. 



7 



to near the threshold region. One possibility then is to choose the transition point (cf. Fig. [I]) at a low value, 
close to m c , so that in effect the 4- flavor calculation by itself covers the full range of physical interest. 

In existing literature, it is already known that the 3-flavor order a 2 s (NLO) calculation can be extended 
to most of the currently accessible energy scales without manifest ill-effects of the large logarithms; and it 
agrees with data rather well. Its efficacy at very large Q is not fully tested (cf. Jll| and results of next 
section). Our calculation will demonstrate the robustness of the complementary, much simpler (order a s ) 
4-flavor calculation. It is worthwhile pointing out that, in the 4-fiavor scheme, the order a s calculation is, 
in fact, also NLO - since the LO 7*c — > c term is of order oP s , as is the case in the standard Q CD t heory of 



inclusive structure functions. This important point is discussed in detail in the Appendix, Sec. A. 5 



3 Inclusive Charm Structure Function 

In the previous section, the theoretical formulation of charm production in DIS is presented within the context 
of the totally inclusive structure functions F\(x,Q). Although terms involving at least one charm quark in 
the final state are the main subject of discussion, they have been added to "light quark contributions" to 
form the totally inclusive structure functions, cf. Eqs. [IT] and |l2|. The reason to do this (rather than simply 
talk about an "inclusive charm structure function" , say F^ ) is: charm quarks are not physically observable; 
there is no unique or obvious definition of a U F£" either experimentally or theoretically. Any experimental 
definition necessarily depends on the procedure or prescription for tagging final state charm (analogous 
to the jet-algorithm for defining jets). Any theoretical definition will be scheme-dependent, and will be 
subjected to questions such as infra-red safety (IRS), e.g. free from ln(m c ) terms. The proper matching of 
the experimental and theoretical definitions is also necessary for a meaningful comparison of theory with 
experiment. 

In this section, we shall be concerned with the comparison of our calculation of inclusive charm production 
structure function with previous calculations, such as that of ACOT and that of the NLO three-flavor 
scheme Q, as well as with existing experimental results. Since the energy range where these comparisons 
will be made is moderate, the question of infra-red safety is not a critical one - as evidenced by the general 
phenomenological success of the NLO 3-flavor calculation (which is not IRS in the sense used in this paper). 
For this purpose, it suffices to define the theoretical 4-flavor F% as the quantity consisting of terms with 
at least one charm quark in the final state in Eq. [ll], i.e. it is the right-hand side less the light-quark 
contribution. At the order we are calculating, this F£ is actually well-defined and infra-red safe, - there are 
no large logarithm terms.| In addition, since the charm mass is just about at the scale for PQCD to become 
valid and, as discussed at the end of the previous section, the transition scale can in practice be chosen close 
to m c , the generalized MS calculation reduces in practice to a ^-flavor (general mass) calculation. 



The 4-flavor NLO calculation: The parton distribution functions /^(x, /i) needed for this calculation 
must be defined in the same renormalization scheme as the hard cross-sections used. In the 4-flavor general 
mass calculation, most CTEQ parton distributions have been defined in the ACOT scheme as described in 
the Appendix, with the matching scales chosen at the heavy quark masses. In particular, the CTEQ5HQ set 
is generated with non-zero heavy quark masses in the hard cross-sections, hence it is the relevant one for our 
application. This set is evolved using NLO evolution, which is appropriate to our calculation which is NLO 
at high energies. For the numerical results presented below, we use an updated version of CTEQ5HQ p^ |. 
In order to compare with results from the LO and NLO 3-flavor scheme, we need the corresponding parton 
distributions in the 3-flavor scheme. The CTEQ5F3 distributions satisfy this need, since they are obtained 
from global analysis of the same data sets as the CTEQ5HQ set, but with parton distributions and hard 
cross-sections in the 3-flavor scheme. 

The calculation of the hard cross-section at 0{a s ) is based on Eqs. j7| and g. In the expression for the 
overall structure functions, the combination j 9 A ® 1 to c g c x + — f\ <8> 1 fjj) ® ^^, due to the subprocesses 

'However, at the next order (ct^) and beyond, F£ is, strictly speaking, not totally infra-red safe by itself (in the sense that 
it does contain some un-cancelled log(^2-) factors); only the sum with the light-quark contributions (in Eqs,p|,[l"l"|) is free from 
such potentially large logarithms. Thus, the 4-flavor formula without the light-parton term should not be used far beyond the 
physical range Q ~ fj, > m c where it is well-defined. This issue is discussed in Ref. []Lq|, where a cutoff parameter is used to 
remove the "soft" parts of the charm contribution from FS and to combine them with the light-parton terms. 
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7*c — » c and 7*5 — > cc, comprise the original ACOT calculation p"2| ]. With non-zero m c , they are all finite. 
The implementation of these terms in the new Monte Carlo calculation is straightforward. We have verified 
that the new MC program reproduces the original ACOT results in detail. Illustrations of the relative 
contributions of the various terms can be found in Ref. [0. Also, the smooth transition of the results of 
the general formalism from the 3-flavor one near threshold to the 4-flavor one at high energies is described 
in that paper. 

The additional contributions f% ® — fx® 1 /<? ® ° w c,A' due to tne T* c ~~ * 9 C subprocess and virtual 

correction to the Born 7*0 — > c term, are new for this calculation. They contain soft divergences which must 
be cancelled. In this Monte Carlo implementation, we use the phase-space splicing method to achieve the 
proper cancellation of the soft divergences between the real and virtual parts. Details of this calculation are 
contained in p5j . For double-checking, we have independently implemented an analytic calculation based 
on the formulas of Hoffmann and Moore ^6|. The two calculations agree quite well with each other over 
the full phase space, with the exception of small values of Q/m c . The difference could be attributed to a 
different treatment of the massive charm quark kinematics adopted by Ref. p6| in deriving their formulas. 
This effect goes away when m c is small compared to Q, as expected]] 

Results and Comparison with order a 2 3-flavor calculation: We now present some typical numer- 
ical results on the theoretical F£(x,Q) obtained in our order a s 4-flavor calculation compared to those of 
order a 2 3-flavor scheme. The 3-flavor results were obtained from the parametrization of Ref. |27|.[] Each 
calculation corresponds to a different way of organizing the perturbation series, hence has its natural region 
of applicability, as discussed in the previous sections. To do a meaningful comparison, it is important to 
take into account the estimated uncertainties of each calculation. Following the example of Fig. |l|, we plot a 
band for each of the predictions, obtained with two common choices of the renormalization and factorization 
scales: /1 = c^J Q 2 + m 2 with c = 0.5 and 2.0. (This choice is more natural than that of the original ACOT 
paper, c.f. footnote |b|.) The 3-flavor calculation requires parton distributions defined in the same scheme. 
We use the CTEQ5F3 set. The results vary in appearance, depending on the kinematic variables. Two rep- 
resentative plots relevant for the HERA measurements are shown in Fig. || where F£(x, Q) is plotted against 
Q for two values of x: (a) x = 0.01; and (b) x — 0.0001. These constitute real examples of the cartoonistic 
Fig. [l], which is designed to emphasize the underlying ideas. These plots show that the overlapping region 





1(T nr «r 
Q»2 (GeV"2) 



Figure 5: Comparison of results of the present order a s generalized MS 4-flavor calculation to that of order a 2 
3-flavor calculations for: (a) x = 0.01; and (b) x = 0.0001. Fg is plotted against Q. The uncertainties of the 
calculation are represented by a band bounded by two curves obtained with the scale choices /1 = cy Q 2 + m 2 ', 
c = 0.5,2. 



■iThe authors of Ref. J2l| have also identified some differences between their recent calculation with Ref. [pr| . We thank S. 
Kretzer for providing us with some details of this comparison. Since the difference in question is numerically not significant, 
none of the results presented below are affected by this problem. 

k We thank Brian Harris for furnishing us with his interface to this program. 



9 



of the two schemes is quite wide for x = 10~ 4 ; and the overall behavior conforms with expectations. For 
x = 1CP 2 , the overlapping region is more limited and it is confined to low Q values. The 3-flavor results fall 
below the 4-flavor ones at large values of Q where the latter should be more reliable. Compared to current 
data, both are within the experimental range; the 4-flavor results are closer to the data points, cf. next 
subsection, although the uncertainties of the 4-flavor calculation are fairly large for x = 10 -4 . 

Since the order a s 4-flavor results are comparable to the order a 2 ones at energy scales close to the 
threshold in both cases (and for other values of x), it is reasonable to choose the transition scale at a 
relatively low value, as mentioned earlier in the paper. The band representing the 3-flavor calculation does 
become wider at large Q for x = 1CP 2 (where the absolute values also are lower than the 4-flavor calculation 
and data); but for x = 1CP 4 , it remains quite narrow. Thus, the theoretically infra-red unsafe logarithms, 
In ' (/i/m c ), do not seem to cause serious problems, at least for very low x. 



Results Compared to recent Zeus data: The general agreement between existing data on charm 
production with order a 2 3-flavor calculations, using the scale choice /i = \J Q 2 + m 2 , is well known. It is 
of interest to compare the same data on "inclusive charm production structure function" F% with our order 
a s 4-flavor calculation. Fig. _| compares the results of our calculation, using the same scale choice, with 
data from the recent ZEUS p8| data. The agreement is clearly excellent. Data also agree with the order a 2 
3-flavor calculation as shown in |2S| ]. This higher order calculation is obviously much more elaborate than 
the order a s 4-flavor calculation presented here. We see that, for inclusive cross-sections, the resummation 
of the a™ln n (/j,/m c ) terms into the charm distribution function f c (x 1 /j,) in the 4-flavor scheme offers a more 
efficient wa y to organize the perturbative series, resulting in an effective NLO calculation already at order 



a s , cf. Sec. A. 5 




Figure 6: Order a s 4-flavor general mass calculation compared to the recent ZEUS data on F$. The 
calculation uses the CTEQ5HQ parton distributions with m c — 1.3GeV. 



4 Semi- inclusive Cross Sections with Tagged Charm Hadrons 

4.1 General considerations 

Now we consider semi-inclusive cross sections, with a charm hadron tagged in the final state. Naively, to 



compute this cross section, one simply convolves the cross sections for parton final states, Eqs. 11 and |12|, 
with a suitable fragmentation function of partons into the final charm hadron. However, the factorization 
of the final state particles through fragmentation functions is only rigorously defined in the limit Q 2 ^> m 2 . 
Thus, the treatment of tagged charm particles in the final state can only be systematically applied at 
high energies, using the 4-flavor scheme. However, it is a common practice to introduce fragmentation 
functions into charm hadrons even in the 3-flavor scheme, and for energies not far above threshold. This 
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approach should be considered a convenient phenomenological model of hadronization, perhaps adequate 
for current experimental accuracy, rather than rigorous theory. We follow this practice in our calculation, 
employing fragmentation functions over the full range of Q 2 , while maintaining the correct factorization- 
scheme implementation at high Q 2 . 

The fragmentation functions d^(x,^) obey the standard (mass-independent) QCD evolution equations, 
and are determined from suitable initial functions at some given scale /x = Qo, of the order of m c . Following 
Mele and Nason , for a given final-state charm hadron H, we write 

df (x, Q ) = d c a (x, Q ) ® D?(x, Q ) (13) 

where the partonic charm fragmentation functions {d c a ] a — l,c} are perturbatively calculable, and D^f(x, Qo) 
is considered non-perturbative and is to be obtained by comparison with experiment. 
For the perturbatively calculable fragmentation functions, Ref. f2i| gives, to order a s : 



d c c (x,Q ) = 6(1 -x) + 



a s (Qo)C, 



F 



2tt 



)2 



1 — x \ m 

Y2 



ln^o _2]n(i_ar)-l 



d c g (x,Q ) = asiQo)TF (x 2 + (l-x) 2 )\n^ (14) 
a z7r 

d q,q,c( X ' Qo) = 

where Tp — 1/2 and Cf = 4/3. In keeping with the choice of the matching scale in our overall calculation, 
we choose Qo = m c for convenience in this paper. 

For the non-perturbative charm quark into charmed mesons fragmentation function (z), we used the 
conventional Peterson form |j(J , 

D " M = .[i-i/Aa-.r <15) 

For the charm meson flW, which will be our focus because of available experimental data, we take e = 0.02, 
cf. plfl , and a value for A such that the branching fraction B(c — » D*) = 0.22 |3l| . We note that, although 
the perturbative fragmentation functions, Eq. [l4|, contain singular (generalized) functions, the overall parton- 
to-charm-meson fragmentation functions d^ (x, fi), Eq. |l3|, are well behaved after convolution with the above 
non-perturbative fragmentation function D® { ) (z). 

In principle, after evolving to high enough Q 2 so that a s \og(Q 2 /m 2 ) is of order one, all of the fragmen- 
tation functions d^ , d^ , d^q 5 eventually become of the same size. In practice, however, at HERA energies 
we find d^ 3> d^ 3> d^q e . For currently required accuracy, it suffices to keep only the charm-to-hadron 
contributions, proportional to d^ . In this approximation, our calculation of the cross section with a tagged 
charm hadron can be written in the 4-flavor scheme as 

Fl x (Q 2 ,x,..) = /S® 1 ^®*? 

+ (I A - Ia® l fg -fl® Vc) ® °<.A ® df (16) 

+ /1® 1 ujI x x ® df . 

To ensure that this calculation is adequate, we have also calculated the contribution from one of the 
more important remaining subprocesses: gluon fragmentation in an order a s light parton hard scattering, 
i.e. j*q — > gq ; g — > H . It is given by: ] q A ® 1 S* 9 A ® dg , where d^ is the gluon fragmentation function 
computed from Eqs. ^ and^. We have verified that its contribution remains small throughout the current 
energy range. It becomes more noticeable only in the large Q limit. However, in this limit, the gluon 
fragmentation function term is not infra-red safe by itself. To insure consistency at high energies, one needs 
to include a full set of infra-red safe higher order subprocesses along with it. The full calculation is more 
appropriately considered as part of the next order project. 

At the same level of accuracy, it is also reasonable to ignore the evolution of d^ (x, fj,(x,Q)), since the 
effect of QCD evolution is not significant over the currently accessible HERA Q range. One can use the 
un-evolved d^(x,/j, = Qo) in place of the fully evolved d^(x, /i(x, Q, m c )) with much gain in efficiency of 
calculation and little sacrifice in accuracy. The error incurred is of the same order as that incurred by 
neglecting the subleading fragmentation functions, dg and d^q S ; and the comments on accounting for 
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ln(m c ) factors made there also apply here. To be sure about this, we have performed the calculation both 
with and without evolving . The difference is indeed small. Therefore, in all our subsequent plots we 
shall include only the direct charm-to-hadron contributions of Eq. (|l^) with the un-evolved fragmentation 
function d^ (x, Qo). 

4.2 Differential distributions 

We employ the Monte Carlo method to carry out the numerical phase-space integration; the new program 
package is implemented in the CH — h programming language. Therefore, we can generate differential distribu- 
tions involving final-state charm mesons, incorporating kinematic cuts appropriate for specific experimental 
measurements, in addition to fully inclusive cross-sections. 

In working with on-mass-shell heavy flavor quarks and hadrons in the parton language, there is an 
ambiguity in defining the momentum fraction variables x (z) for the parton distribution (fragmentation) 
function. This problem arises in all schemes; and it goes away at high energies (Q 3> m c ): where the 
parton picture becomes accurate. Following the modern practice in proofs of factorization, we define the 
momentum fraction variables as ratios of the relevant light-cone momentum components, e.g. p^ = zp£ for 
fragmentation of a charm quark into a D meson. Other authors, e.g. (§[ [j4j, use the prescription pp — zp c 
and adjust the energy variable to enforce the mass-shell condition. At moderate energies, any noticeable 
differences in results due to the choice of this prescription signals that the calculation using fragmentation 
functions is outside the region of applicability of the parton formalism. We have verified that the results 
presented below are insensitive to the choice between the two prescriptions .| 

The QCD formula, Eq. |l^, contains three scale choices in principle: the renormalization scale, the 
factorization scale and the fragmentation scale. For simplicity, we choose the same energy scale fj,(x,Q,m c ) 
for all three. As in the case of the inclusive F^x, Q), for results shown below, we choose the simple function 
/i = c(Q 2 + m 2 ) 1 / 2 , which characterizes the typical virtuality of the process. The constant c is of order 
1; and is varied over same range when we try to estimate the scale-dependence of the physical predictions. 
The magnitude of the charm cross-section is sensitive to the value of m c . For results presented here, we use 
to ms _ \3GeV : the value used in the CTEQ5HQ parton distribution analysis (which is in the middle of the 
range given by the PDG review). 

Fig. shows plots of four differential distributions for D* production at HERA, calculated using the NLO 
(a s ) generalized MS 4-flavor formalism described above. The kinematic variables and their ranges correspond 
to those of the 1996-97 ZEUS data: 1< Q 2 < 600 GeV 2 ; 0.02 < y < 0.7 ; 1.5 < p%' < 15 GeV ; \q D * | < 1.5. 
Each distribution contains two curves obtained with two values of the constant c = 0.5, 1 in the definition 
of the scale parameter described above. These predictions, using the CTEQ5HQ parton distributions, are 
compared to the ZEUS data pj|. We observe a rather large scale dependence in these results. This is not 
surprising, since the compensation among the various subprocesses which underlie the scale-independence 
of the physics predictions (up to some order of perturbation theory), strictly speaking, only apply to the 
inclusive cross-section. The experimental kinematical cuts implemented in these exclusive calculations to 
some extent undermine the mutual cancellation between diagrams which are necessary for relatively scale- 
independent predictions. For example, the order a° j*c — > c HE term (which resums the logarithms arising 
from the near-collinear configurations of an infinite tower of higher-order diagrams) implements the full 
contribution in collinear kinematics, a clear over-simplification. 

Keeping this fact in mind, and with current relatively large experimental errors, the results of Fig. [7] can 
be considered rather encouraging: the Q 2 and pt distributions show very good general agreement; while 
the W and t\d distributions are "in the right ball park", the shapes are too scale-dependent to allow for 
meaningful "predictions" . (A specific choice of scale, in between the two shown, will actually yield theory 
curves in reasonable agreement with data, within errors.) In order to make genuine predictions on differential 
distributions in the 4-flavor scheme, it is necessary to extend the calculation to order a 2 , which would be 
NNLO in the 4-flavor scheme. This can be done by transforming already available NLO results for 3-fiavor 
calculations into the 4-flavor scheme. At the same order in a s , the 4-flavor scheme calculation is, of course, 

'There are some differential distributions, especially those associated with the unphysical partons (such as momentum 
fraction carried by the charm quark, sometimes seen in the literature) which are more sensitive to the choice of definition of 
the momentum fraction variable. 
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more involved than the (NLO) 3-flavor one because of the need for including the necessary subtraction terms 
in a NNLO calculation - such as those appearing in Eqs. ^, ^, and O. 

The calculation of these differential distributions in the a 2 3-flavor FFN scheme was carried out by 
|3^] . Generally good agreement between these calculations and the recent ZEUS data, using specific parton 
distributions, scale choices, etc. has been reported in Ref. [psfl . Although the dependence of the predictions 
to the scale choice was not discussed in this comparison, it is relatively mild, according to [M. This is 
to be expected, because the sensitivity to cuts is reduced with a better approximation of the final state 
particle configurations provided by the order a 2 calculation. This fact implies greater predictive power for 
the differential distributions than the order a s 4-flavor calculation. 

The rjD* distribution in both the 3-flavor and the 4-flavor calculations appear to differ in shape compared 
to the existing data points. This could be due to the inadequacy of applying the fragmentation function 
approach at less than asymptotic region, as discussed in the beginning of this section. In particular, if the 
D* is not collinear to the parton, as assumed in this approach, the rapidity distribution will be affected. 
More extensive study of this effect is obviously needed. 
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Figure 7: Four differential distributions for D* production at HERA, calculated using the NLO (a s ) gener- 
alized MS 4-flavor formalism and m c — 1.3GeV. The two curves on each plot results from two choices of the 
scale constant c 



0.5, 1. The data points are results from the ZEUS 96-97 run 28 



5 Conclusions 

In this work we have shown how the generalized MS formalism can be used to calculate the production of a 
heavy quark over a wide range of energies, from threshold to high Q 2 . For charm production, the formalism 
consists of a 3-flavor scheme calculation at low energies, a 4-flavor scheme calculation at high energies, a 
matching condition between the two schemes, and a transition scale chosen at which one switches between 
the two schemes. Specifically, we have extended the original ACOT calculation for charm production at 
HERA by adding those terms which are necessary to bring the 4-flavor part of the calculation to NLO 



13 



accuracy at high energies. This brings our calculation to the same level of accuracy as the other theoretical 
inputs to the CTEQ and MRS global QCD analyses. 

The generalized MS formalism is ideally suited for inclusive calculations, such as -^2(2;, Q) in DIS. We have 
shown that for a physically-motivated choice of the renormalization/factorization scale, fi, the transition scale 
can be chosen rather low, so that our 4-flavor scheme can be used in the HERA energy range with excellent 
agreement compared to data, and with considerable economy compared to the more elaborate NLO 3-flavor 
calculation. The calculation is much simpler and the resulting program runs faster. Furthermore, the NLO 
corrections are much smaller in the 4-flavor calculation. This suggests that the resummation involved in 
the charm quark distribution function picks up the most important higher-order corrections even at modest 
energies. It also indicates that the perturbation series in the 4-flavor calculation is very well-behaved. 

Our calculation has been implemented as a Monte Carlo program, so that we have also calculated differ- 
ential distributions for exclusive charm final states. By incorporating experimental cuts in the Monte Carlo 
we are able to ensure that our calculation is compared directly with the data, without the need for any 
theoretical extrapolation to all of phase space. The agreement with HERA data is reasonable. However, in 
this case the 3-flavor NLO calculation has somewhat of an edge, because there is no approximation on the 
kinematics as is necessary in the resummation used in the 4-flavor calculation. 

The comparison between available data and the order a s 4-flavor and order a 2 3-flavor calculations 
discussed in the last two sections demonstrates the complementary nature of the two schemes - both with 
regard to the kinematic regions and to the physical quantities they are suitable for. With more abundant 
and more precise experimental data on charm production, the composite generalized MS formalism, which 
encompasses both, will be needed to make reliable comparisons. With this in mind we are now ready 
to extend the 4-flavor calculation to include the necessary 0(a 2 ) terms, along with the 0(a 2 ) matching 
conditions. Such a calculation would include all the advantages the current calculation has in addition to 
the advantages of the current 3-flavor NLO calculation. 

Finally, the generalized MS formalism provides a framework to extract information on the gluon and the 
charm distributions of the nucleon - to the same accuracy as the overall NLO global QCD analysis based 
on total DIS structure functions and other hard processes. Having established that our calculation does a 
reasonable job in describing the existing HERA data, we are now in a position to explore the question of 
whether the proton contains a non-perturbative charm component f35|| . Certainly, it must at some level, 
so the real question is how small is it? The handle on the charm distribution is unique to the generalized 
formalism, since the fixed 3-flavor scheme does not allow the charm parton as an independent degree of 
frcedom.0 We note that there have been recent phenomenological studies of "intrinsic charm" which take 
the theoretical cross-section to be the simple sum of the 3-flavor FFN scheme formulas and an intrinsic charm 



contribution by the heavy-quark excitation mechanism |36, [37|, pq| . This approach cannot be internally 



consistent, because the 3-flavor calculation assumes parton evolution with no charm quark distribution, 
while "intrinsic charm" explicitly requires one. The intimate interplay between the hard matrix elements 
and parton evolution is consistently incorporated only in the generalized MS scheme. Although there has 
been a recent effort to examine this problem incorporating some of the ideas of the general scheme |3!| , a 
fully consistent study, preferably based on more extensive data, still awaits to be done|] 



Acknowledgments 

We would like to thank Brian Harris for many helpful discussions and especially for providing us with the 
code for the NLO 3-flavor calculations, and for verification of the numerical results. We would also like 
to thank John Collins and Fredrick Olness for many discussions and help; Jim Whitmore for assistance 
about the ZEUS data; and Jack Smith for discussions on the extensive order a 2 , calculations carried out by 

m The CTEQ5HQ parton distributions used in our calculations here contains charm partons, but does not use an independent 
non-perturbative charm component as input. This feature is not inherent to the formalism. 

"Because all active partons are coupled, the inclusion of a non-perturbative charm component of the nucleon will affect all 
parton distributions. Hence, it will influence the global fitting of all available data sets, including the extensive DIS sets. If 
a reliable conclusion is to be drawn, it will not suffice to combine a subset of existing- pa rton distributions (in a conventional 
scheme) with modified charm and gluon distributions (in the new scheme), as is done in |39| using GRV and MRST distributions. 
Phcnomenologically, such a combination can upset the original good global fit to existing data, particularly the precision DIS 
data. Among theoretical problems, the NLO evolution equation and the momentum sum rule will not be observed, due to the 
mixing of these hybrid parton distributions. 



14 



the Stony-Brook/Leiden group and for clarification of the proper references to their work. This work was 
supported in part by the US National Science Foundation under grants PHY9802564 and PHY9722144. 



15 



A Formalism for Inclusive DIS Cross Section 



In this appendix, we give a concise and careful presentation of the general formalism, including the relatively 
simple operational procedure for calculating its various components. We hope this will help fill the gap 
between the original, relatively sketchy, ACOT paper (l2) and the recent, more technical, all-order proof of 
the formalism by Collins jl3| . 

The basis for all discussions is the factorization theorem in the presence of non-zero quark masses fl3|j , 
Eq. ||. To establish this theorem in PQCD (cf. Eq. [l?] below), and to give precise meaning to the various 
factors, one must work with partonic cross- sections and parton distributions inside partons, rather than the 
corresponding physical quantities which appear in Eq. In this concise summary, we proceed as follows: (i) 
spell out the operational procedure to establish the factorization formula and calculate the hard cross-sections 
in any scheme; (ii) describe the specifics of the 3-flavor and 4-flavor schemes respectively; (iii) discuss the 
matching conditions between the two; and (iv) consider the transition from one to the other in the composite 
scheme, which constitutes the general formalism of Collins El]. We finish with some remarks on the meaning 
of "LO" and "NLO" calculations in different schemes. 



A.l Procedure to define the factorization scheme and calculate the hard cross- 
sections 

Given the QCD Lagrangian, with non-zero masses for the heavy quarks, one arrives at the general factor- 
ization formula for partonic cross-sections and parton distributions as follows: ^| 

(i) Start with a set of relevant partonic structure functions ui a similar to the left-hand side of Eq. (Q) but 
with on-shell parton targets and calculate them in perturbation theory in a given renormalization scheme 
(i.e. with specific ultra-violet counter-terms) to a given order in a s . The result uj a {Q,x, \,a s (//)) will 
depend on the renormalization scale /i and will contain collinear singularities (represented by =■) as well as 
potentially large logarithm terms of the form (a s ln(- i5 -)) n . (For simplicity we do not differentiate between 
the renormalization scale hr and the factorization scale /if, both of which are taken in practice to be of 
order Q.) 

(ii) Independently, calculate the set of process-independent perturbative partonic distribution functions f% in 
the same renormalization scheme, using either the (moment-space) operator-product expansion or, equiva- 
lently, the (rc-space) bi-local operator definition of the distribution functions. Both ultra-violet and collinear 
singularities appear in this calculation. The ultra-violet singularities are removed by additional counter-terms 
which, along with the coupling constant renormalization counter-terms, define the factorization scheme. The 
result takes the form f^(x, 21a, a s (p)). 

(iii) Confirm that all collinear singularities in the form of ~ terms, appearing in uj a (Q, x, 21a, i, a s (fi)) appear 
in the universal form given in the process-independent functions fa(x, ^jf-, \ ,a s (/i)), so that they can be 
factorized out in the manner of Eq. (||) , 

uj a { — ,x, — , -,a s (/x)) =S~] fa(x, — , -,a s (fx)) <8> £ b ( — ,x, —,a s (fj,)) , (17) 
[i fi e ^ fi e MM 

with u) a being fully infra-red safe in the sense that it is free of all - dependence. In the 4-flavor scheme 
defined below, the functions f^(x, 2la, ^, a s (/i)) will also contain the same large logarithmic terms as 

x, ^la, i, a s (n)), so that these too factorize in Eq. (|Ti| ) with the result that uj a is free of all (a s ln(^-)) n 
terms, and it is well-behaved as m c — > 0. 



(iv) Systematically invert Eq. (17) to solve for the set of finite hard cross-sections uj a , which are then used 
in Eq. (^) for calculating physical structure functions. 

°This procedure may sound familiar because it is based on the basic principles of PQCD. We include it here because the 
details, especially concerning the heavy quark mass dependence, are quite distinct from conventional practices. Thus, it is 
essential to explicitly spell out the steps involved. 
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There are two points to note: (a) The inversion of Eq. ( |l7| ) order-by-order in the perturbation series is 
equivalent to subtracting the singularities contained in from u a ; (b) There is no need to set the quark 
mass(es) to zero anywhere in the above procedure. 

In the following, we apply the above procedure to define the two simple renormalization schemes, involving 
3 or 4 active quark flavors, which underlies the general approach of Refs. J9|, |l3); and combine them to define 
the latter in the subsection under the heading of the generalized MS formalism. These discussions are 
applicable to all orders in perturbation theory. Throughout these discussions, the three quarks {u, d, s} , 
with masses comparable to or less than A, will be referred to as light quarks, and denoted collectively by 
q. The collection of light quarks plus the gluon g will be referred to as light partons, and denoted by I. 
As mentioned earlier, although the formalism applies to all heavy quarks {c, b, t} , we shall use the case of 
charm as a generic representative, for concreteness and clarity - hence the 3- and 4-flavors. Because the real 
charm quark mass m c is not large compared to the on-set of the region of applicability of PQCD, the 4-flavor 
scheme plays a more prominent role in practical applications discussed in the main body of this paper. For 
a heavier quark, the two corresponding schemes and their proper matching, as discussed in the rest of this 
(theoretical) section, will be more relevant. 

A. 2 Three-flavor Scheme 

The 3- flavor scheme is precisely defined by choosing to work with only 3 active quark flavors, consisting of the 
light quarks, and using the subtraction procedure of Ref. Q. The prescription for subtracting ultra-violet 
divergences encountered in the calculation of the partonic structure functions depends on the particle that 
produces the divergence. Broadly speaking, divergences due to the light partons I, are removed using MS 
counter terms, whereas those due to the charm quark c are removed by BPH zero-momentum subtraction 
counter terms. The precise definition can be found in Ref. |l3"| . This ultra-violet subtraction scheme has the 
nice feature that the charm quark explicitly decouples as its mass becomes large. In particular, the operators 
which make up the charm quark distribution function are suppressed by powers of order A 2 /m 2 . Since these 
terms are power-suppressed in the "heavy quark" mass, they are usually excluded from the 3-flavor scheme 
parton picture. 

In practice then the partonic calculations in this scheme are done by considering diagrams where the 
massive charm quark can only appear in the final state, and there are no charm quark distribution functions, 
cf. Fig. @. The light parton distributions always evolve according to the 3-flavor DGLAP equation, for all 
values of the renormalization scale \x — both below and above the heavy quark production threshold. The 
parton distribution functions defined in this scheme will be restricted to the light parton I = {g, q, q} sector, 
and they will be denoted by 3 f l A . In the perturbative calculation, 3 // contains e~ x pole terms which are due 
to collinear singularities. The lowest order (LO, 0(aD) partonic process in which the charm quark appears 
in this scheme is the j*g — » cc "heavy- flavor creation" (HC) process (also known as boson-gluon fusion), 
corresponding to the diagrams of Fig.^a). The associated partonic structure function, denoted by uo c g c , is 
finite. The next-to-leading order (NLO) contribution includes the 1-loop virtual corrections to j*g — > cc (cf. 
Fig.(||b)), plus the real partonic HC processes 7*Z — ► del (cf. Fig.^c)). The collinear divergences which appear 
in the calculation of the 0(a 2 s ) partonic structure functions 3 w°° and 3 wf cl arise from splitting of massless 
light partons in the collinear configuration, and take the form of e^ 1 pole terms, precisely corresponding 
to those appearing in 3 fj mentioned above. That is, the partonic structure functions have the factorized 
structure shown in Eq. (0), and the hard cross-section functions Cj\ will be free from e _1 collinear singularities. 
This 3-flavor scheme is the one used by Ref. to calculate charm production to NLO, i.e. 0(a 2 ).f] 

At high energies the hard cross-sections calculated in this scheme contain powers of \n(Q 2 /m 2 ), as men- 
tioned in the introduction. The perturbative expansion should be accurate at energy scales not too far 
above threshold, or Q 2 ~ m 2 , where \n(Q 2 /m 2 ) is of order 1. However, at high Q 2 ^> m c the perturbative 
expansion parameter is effectively ct s ln(Q 2 /rn 2 ), and the large logarithm factor spoils the convergence of the 
perturbative series. In other words, the "hard cross-sections" uj a defined in this scheme arc finite, but not 
infra-red safe in the limit ^jf- — > 0. 



p To be consistent, the virtual correction to the process 7*g — » qg, which contains a charm quark loop, must also be included 
at this order in the 3-flavor scheme calculation of the total inclusive structure functions. 
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A. 3 Four-flavor scheme with non-zero m c 

In order to better deal with these logarithms at high energies it is more useful to use the 4-flavor scheme, in 
which the renormalization of oj a and is carried out using dimensional regularization and the MS counter 
terms for all Feynman diagrams, while keeping the full quark mass dependence (m c ) of the Lagrangian. 

Charm distribution functions calculated in this scheme, (a = I, c), are not suppressed as in the 3-flavor 
scheme, but contain powers of ln(m c //x), along with possible e _1 poles. Because of the different subtraction 
procedures used in the two schemes, even the light parton distributions 4 // will differ from 3 fj by a finite 
renormalization in general. (We will return to this point later.) Because renormalization constants in the MS 
subtraction procedure are independent of mass, the evolution kernels for the parton distributions will be 
the same as the corresponding ones in the familiar zero-mass 4-flavor case. This is a significant convenience. 
The perturbative parton distribution functions ^f\ have been calculated to NLO in Ref. pc[ |. 

Since charm also has a parton interpretation in this scheme, the set of partonic processes are expanded 
to include those involving charm initial states. The LO partonic process that involves the charm quark in 
the 4-flavor scheme is the 7*c — > c "heavy-quark excitation" (HE) process (Fig.(^a)). NLO charm quark 
contributions in the 4-flavor scheme come from the 1-loop virtual corrections to HE 7*c — > c (Fig.(|Jb)), 
and from the real HE 7*0 — > gc and HC 7*5 — > cc processes (Fig.(|p,d)). Partonic structure functions 
ui a calculated beyond LO in this subtraction scheme contain both e _1 poles (due to collinear singularities 
associated with light degrees of freedom) and powers of mass-logarithms, ln(Q/m c ), (due to collinear con- 
figurations associated with the heavy degree of freedom), just as in the 3-flavor scheme. The important 
difference compared to the latter case is that these potentially large logarithm terms also appear in the 
4-flavor parton distributions 4 f^. Consequently, they are systematically factored out from ui a when we obtain 
the hard cross-sections u> a by inverting the factorization formula Eq. ([l7|) . The charm distribution function 
represents the resummed contribution of all the large (infra-red unsafe) logarithm terms in w„. As a result, 
u) a is free from both types of collinear "singularities" (in quotes since the logarithms become singular only 
in the zero-mass limit). In effect, all logarithm factors ln(Q/m c ) in uj a are replaced by ln(Q/^) in Cj a , (with 
accompanying finite subtractions), and the latter is infra-red safe in the ^ — > limit.^ Thus, the 4-flavor 
scheme has a well-defined high energy limit, and is expected to give a much more reliable description of the 
physics of charm production at large Q than the 3-flavor scheme. 

As formulated above, the hard cross-sections still contain finite charm-mass dependence, i.e. Cj a = 
u} a (x, — , ^jf Being infra-red safe, as m c /Q — > 0, the limit ui a (x, — , ^f-, u)— >a>™ c=0 (a;, — , u) is well de- 
fined. In this limit, the 4-flavor scheme with non-zero charm mass reduces to the conventional zero-mass 
(ZM) 4-flavor parton schemef], as mentioned in the introduction. As emphasized in Ref. [Q, however, the 
factorization of potentially dangerous ln(m c ) terms does not require taking the m c — > limit in the infra- 
red safe coefficient functions. The conventional practice of always setting m = in the hard cross-section 
u> a (x, Q, )Jb) is a convenience, not a necessity; it results from the use of dimensional regularization of the zero- 
mass theory as a simple and efficient way to classify and to remove the collinear singularities. For a "heavy 
quark" with non-zero mass m c , this convenient method of achieving infra-red safety is not a natural one (as 
it is for light flavors), since m c itself already provides a natural cutoff. In other words, the theory has no 
real collinear "singularities" associated with the charm quark, and the universal (i.e. process-independent) 
and potentially large mass-logarithms can be factorized systematically as outlined above. In fact, by keeping 
the charm quark mass dependence, this scheme can be extended down to lower values of Q with much more 
reliable results than in the zero-mass case. This is possible because of the well-defined relation between the 
4-flavor calculation with non-zero m c and the 3-flavor (FFN) calculation; e.g. at order a s , Ref. |l^] showed 
that, for given x, Q 

*FZ( X ,9.,—) V ^ 3 ^-,— ) 0(a s ) (18) 

where the superscripts 3,4 refer to the 3- and 4-flavor scheme calculations respectively. To distinguish this 
more general 4-flavor scheme from the conventional zero-mass (ZM) 4-flavor scheme, we can refer it as the 
general-mass (GM) 4-flavor scheme. 

q The validity of these statements to order ot 2 s can be inferred from the explicit calculations of Refs. j?], 5(1 The proof to 
all orders of perturbation theory has been given in Ref. mM. 

r In conventional zero- mass (ZM) 4-flavor theory, collinear singularities due to charm appear as e~ 1 poles along with those 
from other flavors, and are regulated accordingly. When properly calculated, the massless limit of our (m c 7^ 0) Wilson 
coefficients, u)^ c— (x,Q,fi), should agree with the standard zero-mass results. 
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The theoretical result (Eq. |lj) does not, however, constrain the threshold behavior of the predicted 
physical structure function in the limit of lim Q (orW) — > m c ; to make a physical prediction, one needs 
to first choose ^ as a function of the physical variables {x,Q,m c }. This is related to the well-known scale 
dependence of PQCD prediction in general. We shall return to this problem at the end of the next subsection. 

There is one additional advantage of the 4-flavor scheme. Since the charm quark distribution is explicitly 
included in the 4-flavor scheme, and since m c is not much larger than a typical non-perturbative scale such as 
the nucleon mass, one can allow for the existence of a possible nonperturbative ( "intrinsic" ) charm component 
inside a hadron at a low energy scale, say Qq — as the boundary condition for evolution to higher scales, just 
like the other light flavors. This is a possibility not permitted in the 3-flavor scheme by assumption. 



A. 4 The generalized MS formalism with non-zero m c 

Both the 3-flavor and the 4-flavor schemes described above are valid schemes for defining the perturbative 
series of the inclusive cross section in principle. They are equivalent if both are carried out to all orders in 
the perturbation series. At a given finite order, they differ by a finite renormalization| of the distribution 
functions, as well as the strong coupling a s . From the physics point of view, when calculated to the 
appropriate order (cf. below), the 3-flavor scheme provides a more natural and accurate description of the 
charm production mechanism near the threshold (Q 2 ~ m 2 ), whereas the 4-flavor scheme does the same in 
the high energy regime (Q 2 3> m 2 ), as shown in Fig. |]. 

The precise definitions given in the above subsections provide the means to implement the intuitive ideas 
discussed in the introduction. A unified program to calculate the inclusive structure functions, including 
charm, which maintains uniform accuracy over the full energy range, must be a composite scheme consisting 
of: 

(i) the 3-flavor scheme, applied from low energy scales, of the order of m c , and extended up; 

(ii) the 4-flavor scheme, applied from high energy scales on down; and 

(iii) a set of matching conditions which define the perturbative relation between the two schemes applied 
at a specific matching scale /i m . 

It is useful to explicitly discuss all the elements of this composite scheme which link the component 3-flavor 
and 4-flavor calculations discussed in previous subsections to physics predictions of the general formalism: 
Choice of scale: Within each scheme (i — 3 or 4), one needs to specify /i as a function of the physical 
variables in order to make a physical prediction, i.e. 

F^ s (Q 2 ,x,m c ) = F^ x (Q 2 ,x,m c ,fj,{x,Q,m c )) 

Although there is considerable freedom in choosing fi(x,Q,m c ), two conditions should be met so that 
the prediction can be reliable: (i) pL must be of the order of Q or m c so that PQCD applies, and (ii) 
F^\(Q 2 , x, m c , /i) must be relatively stable with respect to variations of [i for the (x, Q) of interest. This 
is the well-known scale-dependence of any PQCD calculation. In Fig. |l|, the presence of the uncertainty 
associated with the choice of n(x, Q, m c ) in each scheme is represented by the respective bands. 

Matching conditions and choice of matching scale: For a given set of arguments, F^\(Q 2 ,x,m c , /i) 
and F^ x (Q 2 ,x,m c , /i) are not independent. Being the same physical quantity calculated in two different 



schemes (cf. Sec. A.2 and A.3 ), they are related by a finite renormalization: 

4 a s ( M ) = 3 a s {fi) + Aa _ , . 

*f a (x, M ) = 3 f a (x, M ) + Af a (x, M ) PP at, '-' im (19) 

where Aa and Af a (x, /j,) are fully calculable once the two schemes are defined. They have been calculated 
to order a 2 M . Specifically, the simpler results at order a s are H : 



^ In ^ + 0(a 2 s ) 



Gtt 



(20) 



s The magnitude of the "finite" renormalization depends on the renormalization scale: e.g. ln(m c //i) factors are finite, but 
can be numerically large if fi 2> ra c . 
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4 / 9 (^) = 3 /W) + + 0(a 2 s ) 

4 / s (^M) = 3 / 9 (^) ^ln^ 3 / 9 (^M) + 0{al) (21) 

4 / c (^) = + !^ m ^ ff (2 2 + (l-*)V/ 9 (§,M) + 0(a 2) 

The scale at which these two schemes are matched will be called the matching point, and denoted by 
fi m . Note that either scheme can still be used with fi above or below the matching point, it is just that 
the equations ( |l9|) are only enforced at /i m . In principle, /i m can be chosen at any value - different 
choices lead to the same overall results, up to higher order corrections. As can be seen in Eqs. [2(] and 
in the generalized MS scheme, Aa(/Lt) and Af a (x,/j,) are both of the form a s {p) ln(/x/m c ) C\ + 0{a 2 ). 
Thus, if one chooses fj, m = m c , both functions a(n) and f a (x,jj) in the 3-flavor scheme are equal to their 
counterpart in the 4- flavor scheme to first order in a s at the matching point. Most recent works adopt 
this choice || fy| |l8); we do the same in this paper. Although this choice is convenient, it is not required 
in the general formalism. The ideas behind the matching conditions, Eq. |l9|, are illustrated in Fig. || which 
show two possible matching points, the first being the special one fi m = m c . This plot also shows that 
fi m should not be chosen too far above m c , lest the factor a s ln(/i/m c ) in the discontinuity ceases to be 
perturbative. 



OC s (u) or f(x,\i) 

4-flavor 




Figure 8: Schematic illustration of the matching between 3-flavor and 4-flavor schemes. Two possible 
matching points are shown: the first one, /it m i, is at the the value m c where the discontinuity is of order a 2 ; 
the second one, fi m 2, is at an arbitrary value where the discontinuities are of order a s ln(p./m c ). 



Choice of transition scale: With F%\{Q 2 , x, m c , fi(x, Q, m c )), i = 3,4, we have two sets of calculations, 
one for each scheme. For physics applications, wc need to specify which of these to use, say 

pP hy S(n 2 , m s_ r F [ ^ x {Q 2 ,x,m c ,^(x,Q, mc )) Q < Mt , 00 , 
f a,\\Q ,x,m c ,n(x, Q,m c )) Q > \h 

where we have introduced another scale fi t - the transition point - where one switches from one scheme 
to the other, according to which one is more appropriate, as discussed in the introduction, cf. Fig.|l| and 
the previous subsections of this appendix. 

Conceptually, the transition point is distinct from the matching point, as should be clear from their defining 
equations, [l^ and p2|] *| The guiding principle for choosing the optimal /i t is that it should be within the 
region where the and F^ x calculations are both valid, and that their differences within this region 



* This distinction was first made in Collins' paper on the general formalism 
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are small. For instance, in the idealized situation depicted in Fig. |l| \i t is best chosen to be in the middle 
of the Q range where both uncertainty bands are relatively narrow. 

In practice, one can only estimate the range of uncertainties of the 3-flavor and 4-flavor calculations, say 
by examining the scale-dependence of the respective calculations and then make a judicious choice of fi t . 
In the case of inclusive charm production discussed in Sec. ||, Fig. |^ shows that the transition scale can be 
chosen at a relatively low value, close to m c . For this choice, the composite scheme calculation reduces, 
in practice, to just the 4-flavor calculation. 

A. 5 What do "LO" and "NLO" mean? 

As already mentioned in the body of this paper, in a multi-scale problem such as heavy quark production, the 
designation of "LO" and "NLO" to a given calculation can be rather misleading in conventional fixed-order 
calculations, due to the presence of large logarithms which vitiates the naive counting of powers of a s . In 
the composite scheme, which has a wider range of applicability than FFN schemes, the meaning of "LO" 
and "NLO" can be better defined, provided the relative magnitudes of the large scales are properly kept in 
mind. We elaborate a little bit. 

• In the 3-flavor scheme, "LO" consists of the 0{a\) HC j*g — > cc process; whereas "NLO" involves 0(a 2 ) 
processes such as j*g — ► gcc. This formal designation makes physical sense only in the threshold region. 

• In the 4-flavor scheme, the "LO" process is represented by the 0(oP s ) HE 7*0 — » c process; and the "NLO" 
ones consist of the 0(a\) HC process as well as the 0(a\) HE 7*0 — > gc process. This formal designation 
coincides with the familiar one in the conventional treatment of DIS structure functions; it makes physical 
sense only when Q 2 3> m 2 . 

The apparent mismatch of orders in a s (e.g. order a s being "LO" in the former, but "NLO" in the latter) 
can be understood within the composite scheme which takes into account the order of magnitudes of the 
other relevant quantities in the factorization formula at the appropriate energy ranges. 

Specifically, the 0{a\) and 0(a 2 ) 3-flavor calculations lose their "LO" and "NLO" meaning as Q 2 
becomes very large compared to m 2 , since each power of a s is accompanied (and neutralized) by a large 
logarithm ln(Q 2 /m 2 ) factor in the hard cross-section. Consequently, all terms become effectively 0(a? s )\ 
In fact, the resummation of these large logarithms to all orders of a s gives rise precisely to the charm 
distribution, resulting in the 0{oP s ) charm excitation process 7*0 — > c of the 4-flavor scheme. Conversely, in 
the 4-flavor scheme, although the charm quark distribution f^(x,Q) can be considered to be 0(pP s ) at high 
Q 2 (where all flavors are on an equal footing), as we go down to a lower energy range Q 2 ~ m 2 , one finds 
f^(x,Q) ~ a* ln(Q 2 / TO c) ~ C( a s) compared to the dominant gluon and light quark distribution functions 
(assuming there is no large non-perturbative charm component). Therefore, to consistently match the 4- 
flavor calculation onto the LO (O(aD) 3-flavor calculation, one must include both the "LO" 0(pP s ) 7*c — > c 
and the "NLO" 0(a\)j*g — > cc contributions, along with the associated subtraction term in the 4-flavor 
calculation (Cf. Ref. [fijl). This implies that our 4-flavor calculation, which is NLO at high energy scales, 
becomes effectively "LO" near the threshold, because both the formally 0{aP s ) and 0{a\) contributions 
become of the same order of magnitude, 0{a\) - as in the LO 3-flavor scheme. (The calculations in the 
main part of this paper show that, with the natural choice of scale /i = \/Q 2 + m 2 , the numerical predictions 
of this calculation are actually fairly close to those of the NLO 3-flavor calculation, which is consistent with 
this observation because the size of the NLO correction is within the range of uncertainty of a LO result.) 

This mixing of terms with different apparent powers of a s is physically natural ( cf. Fig. ^ ) and logi- 
cally consistent - it is a necessary feature of switching between different primary schemes, since any finite 
renormalization always entails a resummation (i.e. re-organization) of the perturbation series to all orders. 
In more concrete terms, the need for mixing terms of different apparent powers of a s arises when: 

(i) the LO diagrams for different subprocesses start at different orders of a s ; 

(ii) the associated parton densities are of different numerical orders of magnitude (such as between g,q,c); 

(iii) the order of magnitude of a parton distribution changes as it evolves with Q (such as for c in the region 
above the threshold); and 

(iv) the hard cross-section contains logarithms of ratios of energy scales which become large.f] 

"For these reasons, to require a naive uniform counting of powers of a s over a wide range of Q, when a composite scheme 
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The generalized MS formalism, by keeping the physical m c , provides the appropriate scheme to describe the 
underlying physical processes in the different regions encountered in heavy quark production. 
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Charm final states in deep inelastic scattering constitute ~ 25% of the inclusive cross-section at small x as 
measured at HERA. These data can reveal important information on the charm and gluon structure of the nucleon 
if they are interpreted in a consistent perturbative QCD framework which is valid over the entire energy range from 
threshold to the high energy limit. We describe in detail how this can be carried out order-by-order in PQCD in 
the generalized MS formalism of Collins (generally known as the ACOT approach), and demonstrate the inherent 
smooth transition from the 3-flavor to the 4-flavor scheme in a complete order a s calculation, using a Monte Carlo 
implementation of this formalism. This calculation is accurate to the same order as the conventional NLO F% 
calculation in the limit >> 1. It includes the resummed large logarithm contributions of the 3-flavor scheme 
(generally known in this context as the fixed-flavor-number or FFN scheme) to all orders of a s \n(m 2 /Q 2 ). For the 
inclusive structure function, comparison with recent HERA data and the existing FFN calculation reveals that the 
relatively simple order-a^ (NLO) 4-flavor (m c ^ 0) calculation can, in practice, be extended to rather low energy 
scales, yielding good agreement with data over the full measured Q 2 range. The Monte Carlo implementation also 
allows the calculation of differential distributions with relevant kinematic cuts. Comparisons with available HERA 
data show qualitative agreement; however, they also indicate the need to extend the calculation to the next order 
to obtain better description of the differential distributions. 
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1 Introduction 



Recent measurements of charm production in deep inelastic scattering (DIS) at HERA [1, 2] have shown that up 
to 25% of the total cross-section at small-x contains charm in the final state. This is within the expectations of 
perturbative QCD based on conventional parton distributions. We are now in a position to utilize this process to 
study details of the production mechanism of heavy quarks in general, and to extract useful information on the 
charm and gluon structure of the proton in particular. 

Conventional perturbative QCD (PQCD) theory is formulated in terms of zero-mass quark-partons. For pro- 
cesses depending on one hard scale Q, the well-known factorization theorem then provides a straightforward proce- 
dure for order-by-order perturbative calculations, as well as an associated intuitive parton picture interpretation of 
the perturbation series. Heavy quark production presents a challenge in PQCD because the heavy quark mass, mH 
(H = c, b, t), provides an additional hard scale which complicates the situation - it requires a different organization 
of the perturbative series depending on the relative magnitudes of mH and Q. 

The two standard methods for PQCD calculation of heavy quark processes represent two diametrically opposite 
ways of reducing the two-scale problem to a one-scale problem, (i) In the conventional parton model approach used 
in many global QCD analyses of parton distributions [3, 4, 5] and Monte Carlo programs, the zero-mass parton 
approximation is applied to a heavy quark calculation as soon as the typical energy scale d of the physical process 
Q is above the mass threshold mH. This leaves Q as the only apparent hard scale in the problem, (ii) In the heavy 
quark approach which played a dominant role in "NLO calculations" of the production of heavy quarks [6, 7, 8], 
the quark H is always treated as a "heavy" particle and never as a parton. The mass parameter mH is explicitly 
kept along with Q as if they are of the same order, irrespective of their real relative magnitudes. 

The co-existence of these two opposite approaches represents an uneasy dichotomy in the current literature. 
On physical grounds, the zero-mass parton picture of heavy quarks should be applicable at energy scales very 
much larger than the relevant quark mass, mH <C Q, whereas the heavy quark approach (often referred to as the 
fixed-flavor-number (FFN) scheme) should be more appropriate at energy scales comparable to the quark mass 
mH ~ Q. The actual experimental regime often lies in between these two extreme regions, where the validity of 
either approach can be called into question. There is, however, a natural way to incorporate both approaches in a 
unified framework in PQCD which provides a smooth transition between the two. This has been formulated in a 
series of papers over the years, Refs. [9, 10, 11, 12, 13] e , which has now been adopted, in different guises, by most 
recent literature on heavy quark production in PQCD. [15, 16, 17, 18] 

To see the basic ideas behind this unified picture, let us focus explicitly on the production of charm (H = c) in 
deep inelastic scattering. All considerations apply to a generic heavy quark. Consider the PQCD calculation of the 
F2(x,Q) structure function which receives substantial contribution from charm production as mentioned earlier. 
The underlying physical ideas are illustrated graphically in Fig. 1 where the charm contribution to this structure 
function, denoted F2 c (x,Q), is plotted as a function of Q at some fixed value of x. Near threshold Q ~ mc, 



Figure 1: Schematic illustration of the uncertainty bands (hence the relative merits) of the 3-flavor and 4-flavor 
scheme calculations of the charm contribution to the inclusive structure function - and how they can be naturally 
merged into a composite scheme. 

it is natural to consider the charm quark as a heavy particle, and to adopt the 3 active-parton-flavor scheme of 
calculation (the "heavy quark approach" ). As Q becomes large compared to mc, the fixed-order calculation in this 
approach becomes unreliable since the perturbative expansion contains terms of the form as n log™ {mc 2 /Q 2 ) at any 
order n, which ruin the convergence of the series — these terms are not infra-red safe as mc — > or Q — > oo. Thus 
the uncertainty of the 3-flavor calculation grows as Q/mc becomes large. This is illustrated in Fig. 1 as an error 

d We use Q as the generic name for a typical kinematic physical scale. It could be Q, W , or pT, depending on the process. 
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band marked by horizontal hashes which is narrow near threshold but becomes ever wider as Q/mc increases. On 
the other hand, starting from the high energy end (Q ^> mc), the most natural calculational scheme to adopt is the 
conventional 4-flavor scheme with active charm partons. (In this approach, the infra-red unsafe large logarithms 
mentioned earlier are "resummed" and absorbed into the finite charm parton distributions.) However, as we go 
down the energy scale toward the charm production threshold region, the 4-fiavor calculation becomes unreliable 
because the approximation mc = deteriorates as Q — > mc. The uncertainty band of such a calculation is outlined 
in Fig. 1 by the vertical hashes - it is narrow at high energies, but becomes increasingly wider as one approaches 
the threshold region. 

The intuitive ideas embodied in Fig 1 illustrate that: (i) these two conventional approaches are individually 
unsatisfactory over the full energy range, but are mutually complementary; and (ii) the most reliable PQCD pre- 
diction for the physical F2(x, Q), at a given order of calculation, can be obtained by utilizing the most appropriate 
scheme at that energy scale Q, resulting in a composite scheme, as represented by the cross-hashed region in Fig. 
1. The use of a composite scheme consisting of different numbers of flavors in different energy ranges, rather than 
a fixed number of flavors, is familiar in the conventional zero-mass parton picture. The new formalism espoused in 
Refs. [9, 10, 11, 12] provides a quantum field theoretical basis [19, 13] for this intuitive picture in the presence of 
non-zero quark mass. 

This generalization brings about several important distinguishing features and insights. First, after the hard 
cross-section is rendered infra-red safe by factorizing out the In (mc) terms, the remaining finite dependence on mc 
can be kept in the hard cross-sections to maintain better accuracy in the intermediate energy region. f One can then 
show that the unified formalism reproduces the two conventional approaches in well-defined ways in their respective 
regions of applicability [9]-[13]. Secondly, it has become clear recently that there is much inherent flexibility in 
the choice of the transition energy scale (cf. Fig. 1) as well as the detailed matching condition between the 3- and 
4-flavor calculations. This makes it possible to have a variety of different implementations of the new formalism 
[15, 16, 17, 18], with different emphases. Properly understood, this feature can be exploited to lend strength to 
the formalism; but, during the development of the theory, the subtleties have given rise to much misunderstanding 
and confusion in existing literature. 

In this paper, we study charm production in DIS using the general mass formalism, and compare the results 
with recent HERA data. The main goals of this work are: 

(i) Through this specific example, we give a concise and careful presentation of the general formalism, including 
the relatively simple operational procedure for calculating its various components at any order of a s . We hope 
this will help fill the gap between the original, relatively sketchy, ACOT paper [12] and the recent, more technical, 
all-order proof of the formalism by Collins [13]. 

(ii) We carry out the numerical calculations to make concrete the intuitive ideas illustrated in Fig. 1; and to 
demonstrate the validity of the physical principles underlying the composite scheme. 

(iii) We show that the flexibility of the general formalism mentioned above can result in efficient PQCD calculations 
of inclusive quantities at relatively low order in a s compared to FFN calculations - because the relevant physics 
has been effectively "resummed" by the appropriate scheme adopted for the given energy scale. 

Item (i) is basic; it forms the foundation for the other two parts. However, since this part is about the 
clarification of the existing theoretical formalism rather than new work, and since not all readers are concerned 
with theoretical precision, we have elected to place it in the Appendix. Hopefully, this will increase the accessibility 
of the main body of the paper. In regards to the underlying physics, we believe the discussion in the appendix 
should make a useful contribution to clear the confusion and misunderstanding among the various approaches that 
have been proposed in the recent literature following [12]. 

Based on the terminologies discussed in this Appendix, in Sec. 2 we describe the complete order a s calculation 
carried out in this paper in relation to previous work on this subject. The new calculation extends the validity of 
the original ACOT results to NLO in the high energy regime - on the same level as the conventional zero-mass 
total inclusive structure functions. In addition, the new perspective, as illustrated in Fig. 1 above and discussed 
quantitatively in the paper proper, allows a re-assessment of the physical predictions of the order a s calculation 
near the threshold region, making it a viable alternative to the order a 2 . FFN calculation. In Sec. 3 we present 
the numerical results on inclusive charm production, and demonstrate that the validity and the efficiency of the 
general formalism, as described in (ii) and (iii) above, are indeed seen at this order. We show that very good 
agreement with recent HERA data on inclusive F2 C is obtained in practice. For this study, we have developed 
a new implementation of the generalized MS formalism using Monte Carlo methods. This implementation allows 
the computation of differential distributions, with kinematic cuts, such as da D /dpt, da D /dQ, and da D jdr\ 
which we present in Sec. 4. These results, are in qualitative agreement with available data from HERA. However, 
they are not as good as those of the order as 2 calculation in the 3-flavor scheme. This is to be expected for 
differential distributions with experimental cuts, since the (resummed) low-order calculation contains more severe 

f By contrast, in the standard literature the mc — ► limit is tied to the proof of factorization [20] in the first place. This association 
is not needed, cf. [131 



approximations to the kinematics of the final state partons. For future quantitative studies, the general formalism 
needs to be expanded to incorporate higher-order results (adaptable from existing FFN calculations). This point 
is discussed in the concluding section, along with other observations. 

2 Total Inclusive Structure Functions in the general formalism 

We consider the inclusive DIS structure functions, such as F 2 , focusing on the contribution of a massive quark. 
For definiteness, we assume that the only relevant quark with non-zero mass is the charm quark. The generic 
leptoproduction process is depicted in Fig. 2: 

lx+A — > li + X, (1) 

where A is a hadron, £i j2 are leptons, and X represents the summed-over final state hadronic particles. Note that 
X may or may not contain a visible heavy-flavor hadron. After the calculable leptonic part of the cross-section 




Figure 2: The generic inclusive lepto-production process. We are primarily interesting in final states which contain 
some heavy quarks, in particular charm. 

has been factored out, we work with the hadronic process induced by the virtual vector boson 7* of momentum q 
and polarization A: 

>y*(q,\)+A(P)^X(P x ). (2) 

Although our considerations apply to DIS processes induced by W and Z as well, we shall explicitly refer to the 
neutral current interaction with the exchange of a virtual photon 7* in order to be concrete. The cross-section is 
expressed in terms of the hadronic tensor 

W XrT {q, P) = ^ E < P l e A • J i \Px)(2n) i S^ (P + q- P x ){P x \e a ■ J\P). (3) 

X(Px),spin 

where J] denotes a sum over all final hadronic states. In most cases, it suffices to consider the diagonal elements 
of the tensor F x = W xx (q,P). 

The factorization theorem in the presence of non-zero quark masses - assumed in [12] and established to all 
orders in PQCD [13] - states that the inclusive cross-section can be written as a convolution: 

F%(Q\x,..) = J2 fl(x,») ®Z a , x (x,Q/»,Q/m c ,a s (»)) + o(a/qJ (4) 

a 

where f% is the distribution of parton a inside the hadron A, Q a ^ x is the perturbatively calculable hard cross-section 
for 7* -|- a — > X, p is some positive number, /i denotes collectively the renormalization and factorization scales, and 
a convolution in the x variable is implied. The helicity structure functions F x are simply related to the familiar 
Fi,2,3 [11]. 

The exact way that the physical structure function factorizes into the long-distance (f%) and the short-distance 
(oJ a ,\) pieces on the right-hand side of Eq. (4) depends on the scheme used to define the parton distributions. The 
physical structure function F\ should be independent of any calculational scheme; therefore, the definition of the 
hard cross-sections is determined by the subtraction procedure used to define the parton distributions. As discussed 
in the Introduction, the general formalism consists of the 3-flavor scheme at low energy scales, the 4-flavor scheme 
at high energy scales, and a suitably chosen transition region where matching conditions between the two schemes 
are applied. The precise description of these elements of the formalism is given in the Appendix. Operational 



2.1 Previous calculations 



To put the current calculation in context, we first summarize the existing calculations of leptoproduction of charm 
using the precise definitions given in the Appendix. 

• NLO 3-flavor (3a 2 ) calculation [7]: Most dedicated calculations of heavy flavor production in recent years 
have been carried out in this scheme. The LO process is 0(a\) heavy-flavor creation (HC), 7*5 — > cc. The NLO 
processes consist of the 0(a 2 ) virtual corrections to 7**7 — > cc as well as the real HC 7*/ — > ccl process, where I 
denotes any light parton. Cf. Fig. 3. This calculation becomes questionable when Q 3> m c - i.e., it ceases to be 
"NLO" in accuracy, as indicated in Fig. 1, because the perturbative expansion is actually in a s ln(Q 2 /m 2 ,) for 
large Q. 




a b c 



Figure 3: Representative Feynman diagrams contributing to the calculation of the partonic structure functions for 
charm production in the 3-flavor scheme due to: (a) order a] heavy-flavor creation (HC) or gluon fusion mechanism; 
(b) order a 2 virtual correction to (a); and (c) order a 2 real corrections involving gluon and light-quark initial states. 



Zero-mass 4-flavor (ZM) (4a* 



m. 



0) calculation: This is the formalism used in most conventional QCD 



parton model calculations and popular Monte Carlo programs. It represents an approximation to the general 



mass (GM) 4-flavor scheme by setting m c = in the hard cross-sections u a (x, 



n ' Q 



m c /Q->0 



The LO contribution consists of the O(a s ) 7*c — > c heavy-flavor excitation (HE) process. The NLO contribution 
consists of the 0(a\) virtual corrections to 7*c — > c plus the real HC 7*3 — > cc and, 7*c — > gc HE processes. Cf. 
Fig. 4. This calculation is unreliable in the threshold region, as mentioned in the introduction. 





Figure 4: Representative Feynman diagrams contributing to the calculation of the partonic structure functions for 
charm production in the 4-flavor scheme due to: (a) order a° s heavy-flavor excitation (HE) mechanism; (b) order 
a\ virtual correction to (a); (c) order a\ real corrections to the HE mechanism; and (d) order a\ heavy-flavor 
creation mechanism. 

• LO generalized MS calculation (ACOT) [12]: This represents the simplest implementation of the generalized 
formalism [12]. It emphasizes the overlapping physics underlying the 3-flavor and 4-flavor calculations in the 
region not far above threshold. The 4-flavor calculation consists of HE 7*c — > c (Fig. 4a) plus m c ^ HC 
7*3 — > cc (Fig. 4d), with the requisite subtraction term which removes the m c mass-logarithm. Mathematically, 
the result on F 2 (a; ! Q ! /i) can be shown to match that of the 0(a\) 3-flavor scheme calculation (HC 7*3 — > cc) 
as the (unphysical) factorization scale fx approaches m c from above [12]. Physically, the predicted behavior of 
F 2 C for Q ~ m c will depend on the choice of /i as a function of the physical variables, cf. next section. At a high 
energy scale, /x ~ Q ^> m c , this calculation only approximates the 4-flavor NLO results: it contains the most 
important 0(a\) term, HC 7*5 — > cc (because of the large gluon distribution and the need for matching), but 



does not include the smaller 0(a\) terms represented by Fig. 4b, c. This calculation has been further studied by 
Kretzer and Schienbein [21] and Kramer et al [22]. 

• Variations on the "variable flavor number" theme: In recent years, the general approach proposed in 
[9, 12] has been adopted by other groups, starting from different historical perspectives, [15, 16, 17, 18]. In 
contrast to the fixed flavor-number (FFN) approach, these are usually referred to as being in the variable flavor- 
number (VFN) scheme. This terminology has caused some confusion, since although the common theme is that of 
[9, 12], as shown in Fig. 1, the various implementations differ considerably. Whether a particular implementation 
is self-consistent within the general framework of PQCD, or whether two implementations are compatible within 
the accuracy of the perturbative expansion, are often unclear or controversial (e.g. [23]) because of the complexity 
of the multi-scale problem and because of possible misunderstandings. It is beyond the scope of this paper to 
critically review these approaches. By pursuing the goals described in the introduction, we hope to provide a 
clearer picture of the general formalism of [13], hence a better basis to help address some of the controversial 
issues in the future. 



2.2 The full order a s generalized MS calculation 

The calculation of charm contribution to the inclusive structure functions reported in this paper completes the order 
a s calculation in the general formalism, initiated in [12], including all the hard processes of Fig. 4. The additional 
terms, although relatively small numerically at current energies, are required to make the 4-flavor part of this 
calculation truly NLO at high energies, so that it becomes equivalent to the conventional zero-mass NLO theory 
used in most modern analyses of precision DIS data. In the remainder of this section we discuss the theoretical 
issues and uncertainties in the order a s 4-flavor calculation at all energy scales, in order to address the issues 
highlighted at the end of the introduction, Sec. 1. 

At order a\ the 3-flavor component of the composite scheme consists of only the (unsubtracted) 0(a).) HC 
7*<7 — > cc process. The result is standard. Therefore, our main calculation concerns the 4-flavor scheme component. 
At order a] the right-hand side of Eq. (4) consists of three terms 

F i J t{(Q 2 ,x,m c ,ii) = f c A <S>°Z c cX 

+ n ® ^ 

+ light — parton contributions , 

where the superscript (4) indicates that this is a 4-flavor calculation, and the hard-scattering cross-sections l u* x 
for the various subprocesses are calculated from the corresponding partonic cross-sections (without the hat) 

according to the procedures described in the Appendix. A description of each of the terms follows: 
• The leading order 7*c — > c (Fig. 4a) partonic cross-section °lo c c a is infra-red safe, thus 

= °<a • (6) 



The 7*g — > cc (Fig. 4d) partonic cross-section contains a single power of ln/i 2 /m^ which can be factorized into 
the charm distribution function by the subtraction [12] 



^^^a- 1 /;® °<a. (7) 

where the cancelling logarithm with mass-singularity resides in the 0(a s ) perturbative parton distribution func- 
tion 

l fl = ( a .l2ir)P a ^ q {x)\n( l ?lml) (8) 

• The virtual correction to 7*c — > c (Fig. 4b) plus the real 7*c — > gc (Fig. 4c) partonic process also contain 
In (/x 2 /™ 2 ) terms which are factorized into the charm distribution function by the subtraction 

ir-.cx _ 1, ,cX o, ,c 



w c,A - w c,A ~ Jc ^ w c,\ > 



fc ® V,A > (9) 



where the logarithm appears in the 0{a\) perturbative parton distribution function g 

(10) 



Jc ~ 2tt3 



i±^)fm4-i-2in(i-,; 

1 - x / V mi 



g We have calculated these terms keeping a finite charm quark mass. As discussed in Ref. [13, 22], it would also have been consistent 
to calculate diagrams with an initial state charm quark using a zero charm quark mass. The errors near threshold in both methods of 



Substituting Eqs. (6)-(9) into Eq. (5), the right-hand side can be re-organized as 



F^ x (Q 2 ,x,m c ,») = ft® 1 ^ 



+ (fl ~ f 9 A ® X n ~ fl ® 'fc) ® °<a 

+ fl ® 1 



c,A 

light — quark terms 



(11) 



where the ln(/j/m c ) terms in the w ai >, factors are kept intact, and the needed subtraction terms are explicitly 
grouped with the leading 2— >-l term with the same kinematics. This is the form we use for the actual numerical 
calculations, which we implement using a Monte Carlo approach. 

It is useful to compare this calculation with the NLO 3-flavor calculation [7] which can be written, in the same 
notation, as follows 

F^{(Q 2 ,x,m c ,n) = n® 1 ^ 

+ f 9 A® (12) 

+ light — quark terms 

The term common to the two schemes is O (a\) *y*g — > cc, which appears as the first line in both Eq. 11 and Eq. 12. 
For this comparison, one can consider the rest of the terms in these equations as complementary "corrections" to 
the first term. In particular, the last two terms in the 3-flavor formula Eq. 12 are genuine O (a 2 s ) corrections to 
the common term in the threshold region; hence are commonly referred to as NLO. But at high Q 2 ^> m 2 , these 
terms contain large logarithm factors \nQ 2 /m 2 which vitiates the perturbation expansion. They are not NLO in 
this region. 

In the 4-flavor calculation as organized in the form Eq. 11, the last two lines are also effectively O (a 2 ). This 
is because the distribution f A is effectively O (a\) near threshold, and there is a built-in cancellation between the 
leading order 7*c — > c partonic cross-section and the first subtraction term in Eq. 11, as discussed in detail in 
Ref.[12]. Although these O (a 2 s ) "correction" terms to the O (a).) common term (first line) do not contain the 
full NLO O {a 2 s ) corrections at threshold, they do contain all the O (a 2 ) contributions which are enhanced by 
In (Q 2 /m 2 c ) and quickly dominate as Q 2 increases. In fact, these large logarithmic terms have been resummed to 
all orders in perturbation theory via the DGLAP-evolved charm parton distribution contributions in the second 
and third lines. Therefore, in the large Q 2 region the 4-flavor result represents a true NLO calculation. 

In principle, if carried out to all orders in a s , the 3-flavor and 4-flavor calculations in the form of F A ^ X (Q 2 , x, m c , /j), 
« = 3,4, would give exactly the same prediction for all values of the arguments. The dependence on the scheme (i) 
and the scale (p) arises from the truncation of the perturbation series. Therefore, in order to produce a physical 
prediction F A h ^ s (Q 2 ,x,m c ) and to relate the predictions in the two schemes, a number of additional steps must 
be taken. In the presence of a non-zero heavy quark mass, some of these steps are non-obvious; hence they can be 
the source of confusion. An explicit discussion of these elements of the calculation will make clear the flexibility 
as well as the uncertainties inherent in the formalism. This we do in the Appendix, as part of the more precise 
description of the formalism. Here, we mention only two features which are particularly relevant for the subsequent 
discussions of physical predictions. 

First, within each scheme (i = 3 or 4), one needs to specify [i as a function of the physical variables in order to 
make a physical prediction, i.e. 

F A k \ 3 (Q 2 ' x ' m c ) = F { Jl(Q 2 ,x, m c , fi(x, Q,m c )) 

Although there is considerable freedom in choosing /i(x, Q, m c ), two conditions should be met so that the prediction 
can be reliable: (i) /j must be of the order of Q or m c so that PQCD applies, and (ii) F ( A \(Q 2 ,x,m c ,n) must 
be relatively stable with respect to variations of /j for the (x, <5)-range of interest. This is the well-known scale- 
dependence of any PQCD calculation. For the problem at hand, a common choice for fj,(x,Q,m c ) is \jQ 2 + m 2 : 
it represents the typical virtuality of the internal parton lines in the important subprocesses. The presence of 
the uncertainty associated with the choice of n(x,Q,m c ) in each scheme is illustrated in Fig. 1 by the respective 
bands. 11 

Secondly, as shown in Fig. 1, one needs to identify an appropriate scale at which the 3-flavor and the 4-flavor 
scheme predictions are both reasonable and mutually comparable, so that the transition from one to the other in 
the composite scheme can be made smoothly. In the next section we show that, for the inclusive charm production 



h In the original ACOT paper, a different scale function n(x,Q,m c ) was chosen: it was designed to enforce the condition 
F^(Q 2 , x, m c , n(x, Q, m c )) — ¥ (Q 2 , x , m c , m c ) as Q — ¥ m c . In retrospect, this choice is artificial and unnecessary. The gen- 

eral formalism automatically ensures that, F^(Q 2 ,x,m c ,lJ.) — ¥ (Q 2 , x , m c , m c ) as jit — y m c for given (Q 2 ,x) to the order which 
we are working; but it does not place any restriction on the behavior of F phya (Q 2 ,x,m c ) as Q — ¥ m c in any given scheme. As shown 



cross-section, these conditions can be met over a rather large range of Q, extending down to near the threshold 
region. One possibility then is to choose the transition point (cf. Fig. 1) at a low value, close to m c , so that in 
effect the 4-fiavor calculation by itself covers the full range of physical interest. 

In existing literature, it is already known that the 3-flavor order a 2 s (NLO) calculation can be extended to most 
of the currently accessible energy scales without manifest ill-effects of the large logarithms; and it agrees with data 
rather well. Its efficacy at very large Q is not fully tested (cf. [16] and results of next section). Our calculation will 
demonstrate the robustness of the complementary, much simpler (order a s ) 4-flavor calculation. It is worthwhile 
pointing out that, in the 4-flavor scheme, the order a s calculation is, in fact, also NLO - since the LO 7*c — > c 
term is of order a° s , as is the case in the standard QCD theory of inclusive structure functions. This important 
point is discussed in detail in the Appendix, Sec. A. 5. 

3 Inclusive Charm Structure Function 

In the previous section, the theoretical formulation of charm production in DIS is presented within the context of 
the totally inclusive structure functions F\(x,Q). Although terms involving at least one charm quark in the final 
state are the main subject of discussion, they have been added to "light quark contributions" to form the totally 
inclusive structure functions, cf. Eqs. 11 and 12. The reason to do this (rather than simply talk about an "inclusive 
charm structure function", say F%) is: charm quarks are not physically observable; there is no unique or obvious 
definition of a "F 2 C " either experimentally or theoretically. Any experimental definition necessarily depends on 
the procedure or prescription for tagging final state charm (analogous to the jet-algorithm for defining jets). Any 
theoretical definition will be scheme-dependent, and will be subjected to questions such as infra-red safety (IRS), 
e.g. free from ln(m c ) terms. The proper matching of the experimental and theoretical definitions is also necessary 
for a meaningful comparison of theory with experiment. 

In this section, we shall be concerned with the comparison of our calculation of inclusive charm production 
structure function with previous calculations, such as that of ACOT [12] and that of the NLO three-flavor scheme 
[7], as well as with existing experimental results. Since the energy range where these comparisons will be made 
is moderate, the question of infra-red safety is not a critical one - as evidenced by the general phenomenological 
success of the NLO 3-flavor calculation (which is not IRS in the sense used in this paper). For this purpose, it 
suffices to define the theoretical 4-flavor F% as the quantity consisting of terms with at least one charm quark 
in the final state in Eq. 11, i.e. it is the right-hand side less the light-quark contribution. At the order we are 
calculating, this F 2 C is actually well-defined and infra-red safe, - there are no large logarithm terms. 1 In addition, 
since the charm mass is just about at the scale for PQCD to become valid and, as discussed at the end of the 
previous section, the transition scale can in practice be chosen close to m c , the generalized MS calculation reduces 
in practice to a 4-flavor (general mass) calculation. 

The 4-flavor NLO calculation: The parton distribution functions f A (x,n) needed for this calculation must be 
defined in the same renormalization scheme as the hard cross-sections used. In the 4-flavor general mass calculation, 
most CTEQ parton distributions have been defined in the ACOT scheme as described in the Appendix, with the 
matching scales chosen at the heavy quark masses. In particular, the CTEQ5HQ set is generated with non-zero 
heavy quark masses in the hard cross-sections, hence it is the relevant one for our application. This set is evolved 
using NLO evolution, which is appropriate to our calculation which is NLO at high energies. For the numerical 
results presented below, we use an updated version of CTEQ5HQ [24]. In order to compare with results from 
the LO and NLO 3-flavor scheme, we need the corresponding parton distributions in the 3-flavor scheme. The 
CTEQ5F3 distributions satisfy this need, since they are obtained from global analysis of the same data sets as the 
CTEQ5HQ set, but with parton distributions and hard cross-sections in the 3-flavor scheme. 

The calculation of the hard cross-section at 0(a s ) is based on Eqs. 7 and 9. In the expression for the overall 
structure functions, the combination j a A <E> 1 ^> L g c \ + (f% ~ 1 fg) ® 0uj c A' due *° ^ ne subprocesses 7*c — > c and 

7*<? — > cc, comprise the original ACOT calculation [12]. With non-zero m c , they are all finite. The implementation 
of these terms in the new Monte Carlo calculation is straightforward. We have verified that the new MC program 
reproduces the original ACOT results in detail. Illustrations of the relative contributions of the various terms can 
be found in Ref. [12]. Also, the smooth transition of the results of the general formalism from the 3-flavor one near 
threshold to the 4-flavor one at high energies is described in that paper. 

The additional contributions f A ® — f% ® 1 /,? ® c w c c A , due to the 7*c — > gc subprocess and virtual 

correction to the Born 7*c — > c term, are new for this calculation. They contain soft divergences which must be 
cancelled. In this Monte Carlo implementation, we use the phase-space splicing method to achieve the proper 

'However, at the next order (ctl) and beyond, F% is, strictly speaking, not totally infra-red safe by itself (in the sense that it does 
contain some un-cancelled log(^-) factors); only the sum with the light-quark contributions (in Eqs. 5, 11) is free from such potentially 
large logarithms. Thus, the 4-flavor formula without the light-parton term should not be used far beyond the physical range Q ~ n > m c 
where it is well-defined. This issue is discussed in Ref. [18], where a cutoff parameter is used to remove the "soft" parts of the charm 



cancellation of the soft divergences between the real and virtual parts. Details of this calculation are contained in 
[25]. For double-checking, we have independently implemented an analytic calculation based on the formulas of 
Hoffmann and Moore [26]. The two calculations agree quite well with each other over the full phase space, with 
the exception of small values of Q/m c . The difference could be attributed to a different treatment of the massive 
charm quark kinematics adopted by Ref. [26] in deriving their formulas. This effect goes away when m c is small 
compared to Q, as expected^ 

Results and Comparison with order a 2 3-flavor calculation: We now present some typical numerical results 
on the theoretical F£(x, Q) obtained in our order a s 4-flavor calculation compared to those of order a 2 3-flavor 
scheme. The 3-flavor results were obtained from the parametrization of Ref. [27]. k Each calculation corresponds to 
a different way of organizing the perturbation series, hence has its natural region of applicability, as discussed in the 
previous sections. To do a meaningful comparison, it is important to take into account the estimated uncertainties 
of each calculation. Following the example of Fig. 1, we plot a band for each of the predictions, obtained with 
two common choices of the renormalization and factorization scales: n = C\JQ 2 + m 2 with c = 0.5 and 2.0. (This 
choice is more natural than that of the original ACOT paper, c.f. footnote h.) The 3-flavor calculation requires 
parton distributions defined in the same scheme. We use the CTEQ5F3 set. The results vary in appearance, 
depending on the kinematic variables. Two representative plots relevant for the HERA measurements are shown in 
Fig. 5 where F£(x, Q) is plotted against Q for two values of x: (a) x = 0.01; and (b) x = 0.0001. These constitute 
real examples of the cartoonistic Fig. 1, which is designed to emphasize the underlying ideas. These plots show 
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Figure 5: Comparison of results of the present order a s generalized MS 4-flavor calculation to that of order a 2 
3-flavor calculations for: (a) x = 0.01; and (b) x = 0.0001. F£ is plotted against Q. The uncertainties of the 
calculation are represented by a band bounded by two curves obtained with the scale choices /i = CyJ Q 2 + m 2 .; 
c = 0.5,2. 

that the overlapping region of the two schemes is quite wide for x = 10~ 4 ; and the overall behavior conforms with 
expectations. For x = 10~ 2 , the overlapping region is more limited and it is confined to low Q values. The 3-flavor 
results fall below the 4-flavor ones at large values of Q where the latter should be more reliable. Compared to 
current data, both are within the experimental range; the 4-flavor results are closer to the data points, cf. next 
subsection, although the uncertainties of the 4-flavor calculation are fairly large for x = 10~ 4 . 

Since the order a s 4-flavor results are comparable to the order a 2 ones at energy scales close to the threshold 
in both cases (and for other values of x), it is reasonable to choose the transition scale at a relatively low value, 
as mentioned earlier in the paper. The band representing the 3-flavor calculation does become wider at large Q 
for x = 10~ 2 (where the absolute values also are lower than the 4-flavor calculation and data); but for x = 10~ 4 , 
it remains quite narrow. Thus, the theoretically infra-red unsafe logarithms, In ' 2 (ju/m c ), do not seem to cause 
serious problems, at least for very low x. 

Results Compared to recent Zeus data: The general agreement between existing data on charm production 
with order a 2 3-flavor calculations, using the scale choice fi = \JQ 2 + m 2 , is well known. It is of interest to compare 
the same data on "inclusive charm production structure function" F 2 C with our order a s 4-flavor calculation. Fig. 6 

JThe authors of Ref.[21] have also identified some differences between their recent calculation with Ref.[26]. We thank S. Kretzer 
for providing us with some details of this comparison. Since the difference in question is numerically not significant, none of the results 
presented below are affected by this problem. 

k We thank Brian Harris for furnishing us with his interface to this program. 
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compares the results of our calculation, using the same scale choice, with data from the recent ZEUS [28] data. 
The agreement is clearly excellent. Data also agree with the order a 2 3-flavor calculation as shown in [28]. This 
higher order calculation is obviously much more elaborate than the order a s 4-flavor calculation presented here. 
We see that, for inclusive cross-sections, the resummation of the a™ln n (fi/m c ) terms into the charm distribution 
function f c (x,fi) in the 4-flavor scheme offers a more efficient way to organize the perturbative series, resulting in 
an effective NLO calculation already at order a s , cf. Sec. A. 5. 
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Figure 6: Order a s 4-fiavor general mass calculation compared to the recent ZEUS data on F£. The calculation 
uses the CTEQ5HQ parton distributions with m c = 1.3GeU. 



4 Semi-inclusive Cross Sections with Tagged Charm Hadrons 
4.1 General considerations 

Now we consider semi-inclusive cross sections, with a charm hadron tagged in the final state. Naively, to compute 
this cross section, one simply convolves the cross sections for parton final states, Eqs. 11 and 12, with a suitable 
fragmentation function of partons into the final charm hadron. However, the factorization of the final state 
particles through fragmentation functions is only rigorously defined in the limit Q 2 ~3> m 2 . Thus, the treatment 
of tagged charm particles in the final state can only be systematically applied at high energies, using the 4-flavor 
scheme. However, it is a common practice to introduce fragmentation functions into charm hadrons even in 
the 3-flavor scheme, and for energies not far above threshold. This approach should be considered a convenient 
phenomenological model of hadronization, perhaps adequate for current experimental accuracy, rather than rigorous 
theory. We follow this practice in our calculation, employing fragmentation functions over the full range of Q 2 , 
while maintaining the correct factorization-scheme implementation at high Q 2 . 

The fragmentation functions d^(x,ii) obey the standard (mass-independent) QCD evolution equations, and 
are determined from suitable initial functions at some given scale n = Qo, of the order of m c . Following Mele and 
Nason [29], for a given final-state charm hadron H, we write 

d?(x,Q ) = (r a (x,Q )®D?(x,Q ) (13) 

where the partonic charm fragmentation functions {d c a ;a = l,c} are perturbatively calculable, and D^(x,Qo) is 
considered non-perturbative and is to be obtained by comparison with experiment. 

For the perturbatively calculable fragmentation functions, Ref. [29] gives, to order a s : 

d c c (x,Q ) 

d c g (x,Q ) 

d< q,q,c( X ,Qo) 

where Tp = 1/2 and Cf = 4/3. In keeping with the choice of the matching scale in our overall calculation, we 
choose Qo = m c for convenience in this paper. 
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For the non-perturbative charm quark into charmed mesons fragmentation function D^(z), we used the con- 
ventional Peterson form [30], 

D "^ = z[ i-i/ z -e/(i- z) r (15) 

For the charm meson D^*\ which will be our focus because of available experimental data, we take e = 0.02, 
cf. [31], and a value for A such that the branching fraction B(c — > D*) = 0.22 [32]. We note that, although the 
perturbative fragmentation functions, Eq. 14, contain singular (generalized) functions, the overall parton-to-charm- 
meson fragmentation functions d^ (x, /x), Eq. 13, are well behaved after convolution with the above non-perturbative 
fragmentation function (z). 

In principle, after evolving to high enough Q 2 so that a s \og(Q 2 /m 2 ) is of order one, all of the fragmentation 
functions , d^ , d^ s eventually become of the same size. In practice, however, at HERA energies we find 
df ~3> df ~3> d^ q g . For currently required accuracy, it suffices to keep only the charm-to-hadron contributions, 
proportional to d^ . In this approximation, our calculation of the cross section with a tagged charm hadron can be 
written in the 4-flavor scheme as 

F* x (Q 2 ,x,..) = f A ®^ x ®d H 



+ in - a ® i r 9 -n® 7 c c ) ® °< A ® df (i6) 

+ f%® V c * (g) df 



To ensure that this calculation is adequate, we have also calculated the contribution from one of the more 
important remaining subprocesses: gluon fragmentation in an order a s light parton hard scattering, i.e. )'(/ 4 
gq ; g — > H . It is given by: f A <E> 1 uJg\ <E> d^ , where d^ is the gluon fragmentation function computed from Eqs. 13 
and 14. We have verified that its contribution remains small throughout the current energy range. It becomes more 
noticeable only in the large Q limit. However, in this limit, the gluon fragmentation function term is not infra-red 
safe by itself. To insure consistency at high energies, one needs to include a full set of infra-red safe higher order 
subprocesses along with it. The full calculation is more appropriately considered as part of the next order project. 

At the same level of accuracy, it is also reasonable to ignore the evolution of d^ (x, n(x, Q)), since the effect 
of QCD evolution is not significant over the currently accessible HERA Q range. One can use the un-evolved 
df(x,fi = Qo) in place of the fully evolved d^(x, fi(x, Q, m c )) with much gain in efficiency of calculation and 
little sacrifice in accuracy. The error incurred is of the same order as that incurred by neglecting the subleading 
fragmentation functions, d^ and d^ z ; and the comments on accounting for ln(m c ) factors made there also apply 
here. To be sure about this, we have performed the calculation both with and without evolving d^ . The difference is 
indeed small. Therefore, in all our subsequent plots we shall include only the direct charm-to-hadron contributions 
of Eq. (16) with the un-evolved fragmentation function d^(x, Qo). 



4.2 Differential distributions 

We employ the Monte Carlo method to carry out the numerical phase-space integration; the new program package 
is implemented in the C++ programming language. Therefore, we can generate differential distributions involving 
final-state charm mesons, incorporating kinematic cuts appropriate for specific experimental measurements, in 
addition to fully inclusive cross-sections. 

In working with on-mass-shell heavy flavor quarks and hadrons in the parton language, there is an ambiguity in 
defining the momentum fraction variables x (z) for the parton distribution (fragmentation) function. This problem 
arises in all schemes; and it goes away at high energies (Q ~3> m c ), where the parton picture becomes accurate. 
Following the modern practice in proofs of factorization, we define the momentum fraction variables as ratios of 
the relevant light-cone momentum components, e.g. pj, = zpf for fragmentation of a charm quark into a D meson. 
Other authors, e.g. [8, 34], use the prescription pn = zp c and adjust the energy variable to enforce the mass-shell 
condition. At moderate energies, any noticeable differences in results due to the choice of this prescription signals 
that the calculation using fragmentation functions is outside the region of applicability of the parton formalism. 
We have verified that the results presented below are insensitive to the choice between the two prescriptions. 1 

The QCD formula, Eq. 16, contains three scale choices in principle: the renormalization scale, the factorization 
scale and the fragmentation scale. For simplicity, we choose the same energy scale /i(x,Q,m c ) for all three. As in 
the case of the inclusive F%(x, Q), for results shown below, we choose the simple function fi = c(Q 2 + m 2 ) 1 1 2 , which 
characterizes the typical virtuality of the process. The constant c is of order 1; and is varied over same range when 
we try to estimate the scale-dependence of the physical predictions. The magnitude of the charm cross-section is 
sensitive to the value of m c . For results presented here, we use m^ s = 1.3GeV, the value used in the CTEQ5HQ 
parton distribution analysis (which is in the middle of the range given by the PDG review). 



'There are some differential distributions, especially those associated with the unphysical partons (such as momentum fraction 
carried by the charm quark, sometimes seen in the literature) which are more sensitive to the choice of definition of the momentum 



Fig. 7 shows plots of four differential distributions for D* production at HERA, calculated using the NLO (a s ) 
generalized MS 4-flavor formalism described above. The kinematic variables and their ranges correspond to those 
of the 1996-97 ZEUS data: 1 < Q 2 < 600 GeV 2 ; 0.02 < y < 0.7 ; 1.5 < < 15 GeV ; \r) D , \ < 1.5. Each 
distribution contains two curves obtained with two values of the constant c = 0.5, 1 in the definition of the scale 
parameter described above. These predictions, using the CTEQ5HQ parton distributions, are compared to the 
ZEUS data [28]. We observe a rather large scale dependence in these results. This is not surprising, since the 
compensation among the various subprocesses which underlie the scale-independence of the physics predictions (up 
to some order of perturbation theory), strictly speaking, only apply to the inclusive cross-section. The experimental 
kinematical cuts implemented in these exclusive calculations to some extent undermine the mutual cancellation 
between diagrams which are necessary for relatively scale-independent predictions. For example, the order a° 
7*c — > c HE term (which resums the logarithms arising from the near-collinear configurations of an infinite tower 
of higher-order diagrams) implements the full contribution in collinear kinematics, a clear over-simplification. 

Keeping this fact in mind, and with current relatively large experimental errors, the results of Fig. 7 can be 
considered rather encouraging: the Q 2 and pr distributions show very good general agreement; while the W and tjd 
distributions are "in the right ball park" , the shapes are too scale-dependent to allow for meaningful "predictions" . 
(A specific choice of scale, in between the two shown, will actually yield theory curves in reasonable agreement 
with data, within errors.) In order to make genuine predictions on differential distributions in the 4-flavor scheme, 
it is necessary to extend the calculation to order a 2 , which would be NNLO in the 4-flavor scheme. This can be 
done by transforming already available NLO results for 3-flavor calculations into the 4-flavor scheme. At the same 
order in a s , the 4-flavor scheme calculation is, of course, more involved than the (NLO) 3-flavor one because of the 
need for including the necessary subtraction terms in a NNLO calculation - such as those appearing in Eqs. 7, 9, 
and 11. 

The calculation of these differential distributions in the a 2 3-flavor FFN scheme was carried out by [34]. 
Generally good agreement between these calculations and the recent ZEUS data, using specific parton distributions, 
scale choices, etc. has been reported in Ref. [28]. Although the dependence of the predictions to the scale choice 
was not discussed in this comparison, it is relatively mild, according to [34]. This is to be expected, because the 
sensitivity to cuts is reduced with a better approximation of the final state particle configurations provided by the 
order a 2 calculation. This fact implies greater predictive power for the differential distributions than the order a s 
4-flavor calculation. 

The rjD* distribution in both the 3-flavor and the 4-flavor calculations appear to differ in shape compared to 
the existing data points. This could be due to the inadequacy of applying the fragmentation function approach at 
less than asymptotic region, as discussed in the beginning of this section. In particular, if the D* is not collinear 
to the parton, as assumed in this approach, the rapidity distribution will be affected. More extensive study of this 
effect is obviously needed. 

5 Conclusions 

In this work we have shown how the generalized MS formalism can be used to calculate the production of a heavy 
quark over a wide range of energies, from threshold to high Q 2 . For charm production, the formalism consists of 
a 3-flavor scheme calculation at low energies, a 4-flavor scheme calculation at high energies, a matching condition 
between the two schemes, and a transition scale chosen at which one switches between the two schemes. Specifically, 
we have extended the original ACOT calculation for charm production at HERA by adding those terms which are 
necessary to bring the 4-flavor part of the calculation to NLO accuracy at high energies. This brings our calculation 
to the same level of accuracy as the other theoretical inputs to the CTEQ and MRS global QCD analyses. 

The generalized MS formalism is ideally suited for inclusive calculations, such as F%(x, Q) in DIS. We have 
shown that for a physically-motivated choice of the renormalization/factorization scale, /i, the transition scale can 
be chosen rather low, so that our 4-flavor scheme can be used in the HERA energy range with excellent agreement 
compared to data, and with considerable economy compared to the more elaborate NLO 3-flavor calculation. The 
calculation is much simpler and the resulting program runs faster. Furthermore, the NLO corrections are much 
smaller in the 4-flavor calculation. This suggests that the resummation involved in the charm quark distribution 
function picks up the most important higher-order corrections even at modest energies. It also indicates that the 
perturbation series in the 4-flavor calculation is very well-behaved. 

Our calculation has been implemented as a Monte Carlo program, so that we have also calculated differential 
distributions for exclusive charm final states. By incorporating experimental cuts in the Monte Carlo we are able to 
ensure that our calculation is compared directly with the data, without the need for any theoretical extrapolation 
to all of phase space. The agreement with HERA data is reasonable. However, in this case the 3-flavor NLO 
calculation has somewhat of an edge, because there is no approximation on the kinematics as is necessary in the 
resummation used in the 4-flavor calculation. 

The comparison between available data and the order a s 4-flavor and order a 2 3-flavor calculations discussed 




Figure 7: Four differential distributions for D* production at HERA, calculated using the NLO (a s ) generalized MS 
4-fiavor formalism and m c = 1.3GeV. The two curves on each plot results from two choices of the scale constant 
c = 0.5, 1. The data points are results from the ZEUS 96-97 run [28]. 



kinematic regions and to the physical quantities they are suitable for. With more abundant and more precise 
experimental data on charm production, the composite generalized MS formalism, which encompasses both, will 
be needed to make reliable comparisons. With this in mind we are now ready to extend the 4-flavor calculation to 
include the necessary 0(a 2 s ) terms, along with the 0(a 2 s ) matching conditions. Such a calculation would include all 
the advantages the current calculation has in addition to the advantages of the current 3-flavor NLO calculation. 

Finally, the generalized MS formalism provides a framework to extract information on the gluon and the charm 
distributions of the nucleon - to the same accuracy as the overall NLO global QCD analysis based on total DIS 
structure functions and other hard processes. Having established that our calculation does a reasonable job in 
describing the existing HERA data, we are now in a position to explore the question of whether the proton 
contains a non-perturbative charm component [35]. Certainly, it must at some level, so the real question is how 
small is it? The handle on the charm distribution is unique to the generalized formalism, since the fixed 3-flavor 
scheme does not allow the charm parton as an independent degree of freedom." 1 We note that there have been 
recent phenomenological studies of "intrinsic charm" which take the theoretical cross-section to be the simple sum 
of the 3-flavor FFN scheme formulas and an intrinsic charm contribution by the heavy-quark excitation mechanism 
[36, 37, 38]. This approach cannot be internally consistent, because the 3-flavor calculation assumes parton evolution 
with no charm quark distribution, while "intrinsic charm" explicitly requires one. The intimate interplay between 
the hard matrix elements and parton evolution is consistently incorporated only in the generalized MS scheme. 
Although there has been a recent effort to examine this problem incorporating some of the ideas of the general 
scheme [39], a fully consistent study, preferably based on more extensive data, still awaits to be done." 

m The CTEQ5HQ parton distributions used in our calculations here contains charm partons, but does not use an independent 
non-perturbative charm component as input. This feature is not inherent to the formalism. 

"Because all active partons are coupled, the inclusion of a non-perturbative charm component of the nucleon will affect all parton 
distributions. Hence, it will influence the global fitting of all available data sets, including the extensive DIS sets. If a reliable conclusion 
is to be drawn, it will not suffice to combine a subset of existing parton distributions (in a conventional scheme) with modified charm 
and gluon distributions (in the new scheme), as is done in [39] using GRV and MRST distributions. Phenomenologically, such a 
combination can upset the original good global fit to existing data, particularly the precision DIS data. Among theoretical problems, 
the NLO evolution equation and the momentum sum rule will not be observed, due to the mixing of these hybrid parton distributions. 
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A Formalism for Inclusive DIS Cross Section 



In this appendix, we give a concise and careful presentation of the general formalism, including the relatively 
simple operational procedure for calculating its various components. We hope this will help fill the gap between 
the original, relatively sketchy, ACOT paper [12] and the recent, more technical, all-order proof of the formalism 
by Collins [13]. 

The basis for all discussions is the factorization theorem in the presence of non-zero quark masses [13], Eq. 4. 
To establish this theorem in PQCD (cf. Eq. 17 below), and to give precise meaning to the various factors, one must 
work with partonic cross-sections and parton distributions inside partons, rather than the corresponding physical 
quantities which appear in Eq. 4. In this concise summary, we proceed as follows: (i) spell out the operational 
procedure to establish the factorization formula and calculate the hard cross-sections in any scheme; (ii) describe 
the specifics of the 3-flavor and 4-flavor schemes respectively; (iii) discuss the matching conditions between the 
two; and (iv) consider the transition from one to the other in the composite scheme, which constitutes the general 
formalism of Collins [13]. We finish with some remarks on the meaning of "LO" and "NLO" calculations in different 
schemes. 

A.l Procedure to define the factorization scheme and calculate the hard cross- 
sections 

Given the QCD Lagrangian, with non-zero masses for the heavy quarks, one arrives at the general factorization 
formula for partonic cross-sections and parton distributions as follows: 

(i) Start with a set of relevant partonic structure functions uj a similar to the left-hand side of Eq. (4) but with 
on-shell parton targets and calculate them in perturbation theory in a given renormalization scheme (i.e. with 
specific ultra-violet counter-terms) to a given order in a s . The result w a (^jj, x, ^jf-, \, ot s (p)) will depend on the 
renormalization scale [i and will contain collinear singularities (represented by ^) as well as potentially large 
logarithm terms of the form (a s ln(-^-)) n . (For simplicity we do not differentiate between the renormalization scale 
[ir and the factorization scale /if, both of which are taken in practice to be of order Q.) 

(ii) Independently, calculate the set of process-independent perturbative partonic distribution functions /* in the 
same renormalization scheme, using either the (moment-space) operator-product expansion or, equivalently, the 
(:E-space) bi-local operator definition of the distribution functions. Both ultra-violet and collinear singularities 
appear in this calculation. The ultra-violet singularities are removed by additional counter-terms which, along 
with the coupling constant renormalization counter-terms, define the factorization scheme. The result takes the 
form ^,±, a. (»)). 

(iii) Confirm that all collinear singularities in the form of \ terms, appearing in w a (^,x, Zjf-, \,a s (//)) appear in 
the universal form given in the process-independent functions f%(x, \,a s (/i)), so that they can be factorized 
out in the manner of Eq. (4), 

w a ( — ,x, — , -,a s (//)) =^2fa(x,— ,-,a s (p))<S>Q b ( — ,x, — ,a s (p)) , (17) 

with u) a being fully infra-red safe in the sense that it is free of all ^ dependence. In the 4-flavor scheme defined below, 
the functions f^(x, \,a s (p)) will also contain the same large logarithmic terms as to a (Q,x, \,a s (p)), so 
that these too factorize in Eq. (17) with the result that Q a is free of all (a s ln(^-))™ terms, and it is well-behaved 
as m c — > 0. 

(iv) Systematically invert Eq. (17) to solve for the set of finite hard cross-sections Q a , which are then used in Eq. (4) 
for calculating physical structure functions. 

There are two points to note: (a) The inversion of Eq. (17) order-by-order in the perturbation series is equivalent 
to subtracting the singularities contained in from uj a ; (b) There is no need to set the quark mass(es) to zero 
anywhere in the above procedure. 

In the following, we apply the above procedure to define the two simple renormalization schemes, involving 3 
or 4 active quark flavors, which underlies the general approach of Refs. [9, 13]; and combine them to define the 
latter in the subsection under the heading of the generalized MS formalism. These discussions are applicable to all 
orders in perturbation theory. Throughout these discussions, the three quarks {u, d, s} , with masses comparable to 
or less than A, will be referred to as light quarks, and denoted collectively by q. The collection of light quarks plus 
the gluon g will be referred to as light partons, and denoted by I. As mentioned earlier, although the formalism 

°This procedure may sound familiar because it is based on the basic principles of PQCD. We include it here because the details, 
especially concerning the heavy quark mass dependence, are quite distinct from conventional practices. Thus, it is essential to explicitly 
SDell out the steDS involved. 



applies to all heavy quarks {c, b, t} , we shall use the case of charm as a generic representative, for concreteness and 
clarity - hence the 3- and 4-flavors. Because the real charm quark mass m c is not large compared to the on-set 
of the region of applicability of PQCD, the 4-flavor scheme plays a more prominent role in practical applications 
discussed in the main body of this paper. For a heavier quark, the two corresponding schemes and their proper 
matching, as discussed in the rest of this (theoretical) section, will be more relevant. 

A. 2 Three-flavor Scheme 

The 3-flavor scheme is precisely defined by choosing to work with only 3 active quark flavors, consisting of the light 
quarks, and using the subtraction procedure of Ref. [19]. The prescription for subtracting ultra-violet divergences 
encountered in the calculation of the partonic structure functions depends on the particle that produces the 
divergence. Broadly speaking, divergences due to the light partons I, are removed using MS counter terms, whereas 
those due to the charm quark c are removed by BPH zero-momentum subtraction counter terms. The precise 
definition can be found in Ref. [13]. This ultra-violet subtraction scheme has the nice feature that the charm 
quark explicitly decouples as its mass becomes large. In particular, the operators which make up the charm quark 
distribution function are suppressed by powers of order A 2 /to;?. Since these terms are power-suppressed in the 
"heavy quark" mass, they are usually excluded from the 3-flavor scheme parton picture. 

In practice then the partonic calculations in this scheme are done by considering diagrams where the massive 
charm quark can only appear in the final state, and there are no charm quark distribution functions, cf. Fig. 3. 
The light parton distributions always evolve according to the 3-flavor DGLAP equation, for all values of the 
renormalization scale /i — both below and above the heavy quark production threshold. The parton distribution 
functions defined in this scheme will be restricted to the light parton I = {g,q,q} sector, and they will be denoted 
by 3 Ja- I n the perturbative calculation, 3 // contains e -1 pole terms which are due to collinear singularities. The 
lowest order (LO, 0{a\)) partonic process in which the charm quark appears in this scheme is the 7*<? — > cc "heavy- 
flavor creation" (HC) process (also known as boson-gluon fusion), corresponding to the diagrams of Fig. (3a). The 
associated partonic structure function, denoted by u c g c , is finite. The next-to-leading order (NLO) contribution 
includes the 1-loop virtual corrections to j*g — > cc (cf. Fig. (3b)), plus the real partonic HC processes j*l — > ccl 
(cf. Fig. (3c)). The collinear divergences which appear in the calculation of the 0(a 2 s ) partonic structure functions 
3 Ug C and 3 iuf cl arise from splitting of massless light partons in the collinear configuration, and take the form of e -1 
pole terms, precisely corresponding to those appearing in 3 ff mentioned above. That is, the partonic structure 
functions have the factorized structure shown in Eq. (4), and the hard cross-section functions a)/ will be free from 
e _1 collinear singularities. This 3-flavor scheme is the one used by Ref. [7] to calculate charm production to NLO, 
i.e. 0(a 2 )v 

At high energies the hard cross-sections calculated in this scheme contain powers of \n(Q 2 /m 2 ), as mentioned 
in the introduction. The perturbative expansion should be accurate at energy scales not too far above threshold, 
or Q 2 ~ to 2 , where ln(<2 2 /TO 2 ) is of order 1. However, at high Q 2 ~3> m c the perturbative expansion parameter is 
effectively a s ln((2 2 /TO 2 ), and the large logarithm factor spoils the convergence of the perturbative series. In other 
words, the "hard cross-sections" Cj a defined in this scheme are finite, but not infra-red safe in the limit ^ — > 0. 

A. 3 Four-flavor scheme with non-zero m c 

In order to better deal with these logarithms at high energies it is more useful to use the 4-flavor scheme, in which 
the renormalization of to a and /* is carried out using dimensional regularization and the MS counter terms for all 
Feynman diagrams, while keeping the full quark mass dependence (m c ) of the Lagrangian. 

Charm distribution functions calculated in this scheme, 4 /„ (a = l,c), are not suppressed as in the 3-flavor 
scheme, but contain powers of ln(m c //t), along with possible e _1 poles. Because of the different subtraction 
procedures used in the two schemes, even the light parton distributions 4 /' will differ from 3 // by a finite 
renormalization in general. (We will return to this point later.) Because renormalization constants in the MS 
subtraction procedure are independent of mass, the evolution kernels for the A f h a parton distributions will be the 
same as the corresponding ones in the familiar zero-mass 4-flavor case. This is a significant convenience. The 
perturbative parton distribution functions have been calculated to NLO in Ref. [40]. 

Since charm also has a parton interpretation in this scheme, the set of partonic processes are expanded to 
include those involving charm initial states. The LO partonic process that involves the charm quark in the 4-flavor 
scheme is the 7*c — > c "heavy-quark excitation" (HE) process (Fig. (4a)). NLO charm quark contributions in the 
4-flavor scheme come from the 1-loop virtual corrections to HE 7*c — > c (Fig. (4b)), and from the real HE 7*c — > gc 
and HC 7*3 — > cc processes (Fig. (4c, d)). Partonic structure functions u a calculated beyond LO in this subtraction 
scheme contain both e _1 poles (due to collinear singularities associated with light degrees of freedom) and powers 
of mass- logarithms, \n(Q/m c ), (due to collinear configurations associated with the heavy degree of freedom), just 



p To be consistent, the virtual correction to the process 7*5 — >- qg, which contains a charm quark loop, must also be included at this 
order in the 3-flavor scheme calculation of the total inclusive structure functions. [401. 



as in the 3-flavor scheme. The important difference compared to the latter case is that these potentially large 
logarithm terms also appear in the 4-flavor parton distributions 4 f%. Consequently, they are systematically factored 
out from to a when we obtain the hard cross-sections u a by inverting the factorization formula Eq. (17). The charm 
distribution function represents the resummed contribution of all the large (infra-red unsafe) logarithm terms in u a . 
As a result, Co a is free from both types of collinear "singularities" (in quotes since the logarithms become singular 
only in the zero-mass limit). In effect, all logarithm factors \n(Q/m c ) in co a are replaced by ln(<3//x) in il) a , (with 
accompanying finite subtractions), and the latter is infra-red safe in the ^ — > limit. q Thus, the 4-flavor scheme 
has a well-defined high energy limit, and is expected to give a much more reliable description of the physics of 
charm production at large Q than the 3-flavor scheme. 

As formulated above, the hard cross-sections still contain finite charm-mass dependence, i.e. u a = u) a {x, 77, ^Pf , /•*)■ 
Being infra-red safe, as m c /Q — > 0, the limit Lo a (x, ^, 2 ^ L , /«)— >-a)™ c=0 (:E, ^, /i) is well defined. In this limit, the 
4-flavor scheme with non-zero charm mass reduces to the conventional zero-mass (ZM) 4-flavor parton scheme 1 , 
as mentioned in the introduction. As emphasized in Ref. [12], however, the factorization of potentially dangerous 
ln(m c ) terms does not require taking the m c — > limit in the infra-red safe coefficient functions. The conventional 
practice of always setting m = in the hard cross-section ili a (x,Q,n) is a convenience, not a necessity; it results 
from the use of dimensional regularization of the zero-mass theory as a simple and efficient way to classify and to 
remove the collinear singularities. For a "heavy quark" with non-zero mass m c , this convenient method of achieving 
infra-red safety is not a natural one (as it is for light flavors) , since m c itself already provides a natural cutoff. In 
other words, the theory has no real collinear "singularities" associated with the charm quark, and the universal (i.e. 
process-independent) and potentially large mass-logarithms can be factorized systematically as outlined above. In 
fact, by keeping the charm quark mass dependence, this scheme can be extended down to lower values of Q with 
much more reliable results than in the zero-mass case. This is possible because of the well-defined relation between 
the 4-flavor calculation with non-zero m c and the 3-flavor (FFN) calculation; e.g. at order a s , Ref. [12] showed 
that, for given x, Q 



where the superscripts 3,4 refer to the 3- and 4-flavor scheme calculations respectively. To distinguish this more 
general 4-flavor scheme from the conventional zero-mass (ZM) 4-flavor scheme, we can refer it as the general-mass 
(GM) 4-flavor scheme. 

The theoretical result (Eq. 18) does not, however, constrain the threshold behavior of the predicted physical 
structure function in the limit of lim Q (ovW) — > m c ; to make a physical prediction, one needs to first choose [i 
as a function of the physical variables {x,Q,m c }. This is related to the well-known scale dependence of PQCD 
prediction in general. We shall return to this problem at the end of the next subsection. 

There is one additional advantage of the 4-flavor scheme. Since the charm quark distribution is explicitly 
included in the 4-flavor scheme, and since m c is not much larger than a typical non-perturbative scale such as the 
nucleon mass, one can allow for the existence of a possible nonperturbative ("intrinsic") charm component inside a 
hadron at a low energy scale, say Qo — as the boundary condition for evolution to higher scales, just like the other 
light flavors. This is a possibility not permitted in the 3-flavor scheme by assumption. 

A. 4 The generalized MS formalism with non-zero m c 

Both the 3-flavor and the 4-flavor schemes described above are valid schemes for defining the perturbative series of 
the inclusive cross section in principle. They are equivalent if both are carried out to all orders in the perturbation 
series. At a given finite order, they differ by a finite renormalization s of the distribution functions, as well as the 
strong coupling a s . From the physics point of view, when calculated to the appropriate order (cf. below), the 
3-flavor scheme provides a more natural and accurate description of the charm production mechanism near the 
threshold (Q 2 ~ m 2 ), whereas the 4-flavor scheme does the same in the high energy regime (Q 2 ^> m 2 ), as shown 
in Fig. 1. 

The precise definitions given in the above subsections provide the means to implement the intuitive ideas 
discussed in the introduction. A unified program to calculate the inclusive structure functions, including charm, 
which maintains uniform accuracy over the full energy range, must be a composite scheme consisting of: 

(i) the 3-flavor scheme, applied from low energy scales, of the order of m c , and extended up; 

(ii) the 4-flavor scheme, applied from high energy scales on down; and 

q The validity of these statements to order oil can be inferred from the explicit calculations of Refs. [7, 40, 18]. The proof to all 
orders of perturbation theory has been given in Ref. [13]. 

r In conventional zero-mass (ZM) 4-flavor theory, collinear singularities due to charm appear as e -1 poles along with those from other 
flavors, and are regulated accordingly. When properly calculated, the massless limit of our (m c 7^ 0) Wilson coefficients, d)™ c ~ (x, Q, /x), 
should agree with the standard zero-mass results. 

s The magnitude of the "finite" renormalization depends on the renormalization scale: e.g. ln(m c /^t) factors are finite, but can be 
numerically large if ji m c ■ 
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(18) 



(iii) a set of matching conditions which define the perturbative relation between the two schemes applied at a 
specific matching scale /j TO . 

It is useful to explicitly discuss all the elements of this composite scheme which link the component 3-flavor and 
4-flavor calculations discussed in previous subsections to physics predictions of the general formalism: 

• Choice of scale: Within each scheme (i = 3 or 4), one needs to specify /j as a function of the physical variables 
in order to make a physical prediction, i.e. 

F%»'(Q 2 ,x, m c ) = F ( ;\(Q 2 ,x, m c , »(x, Q,m c )) 

Although there is considerable freedom in choosing [i{x,Q,m c ), two conditions should be met so that the pre- 
diction can be reliable: (i) /i must be of the order of Q or m c so that PQCD applies, and (ii) F^\(Q 2 , x, m c , /i) 
must be relatively stable with respect to variations of /j for the (x, Q) of interest. This is the well-known scale- 
dependence of any PQCD calculation. In Fig. 1, the presence of the uncertainty associated with the choice of 
fi(x,Q,m c ) in each scheme is represented by the respective bands. 

• Matching conditions and choice of matching scale: For a given set of arguments, F^ X (Q 2 ,x,m c , /i) and 

F^\{Q 2 , x, m c , /j) are not independent. Being the same physical quantity calculated in two different schemes 
(cf. Sec. A. 2 and A. 3), they are related by a finite renormalization: 

4 a s (ii) = 3 a s (ii) + Aa _ 
4 f*( Xj = 3 / a (a . j ^ + A/ a (;Ej /i) applied at /i - /i m (19) 

where Aa and Af a (x,n) are fully calculable once the two schemes are defined. They have been calculated to 
order a 2 [40]. Specifically, the simpler results at order a s are [9]: 



,2 

A a s (p) = 3 a s (fi) 



6tt ml 



(20) 



4 P(x,ti) = 3 P(x, M ) + + <D(a 2 ) 

4 f 9 (x,») = 3 f°(x,n) - S -^\n^ 3 p(x,(i) + 0(a 2 ) (21) 

4 nx,u) = + !^Mln^/f (z 2 + (l-z) 2 ) 3 / fl (f,^) + 0(a 2 ) 

The scale at which these two schemes are matched will be called the matching point, and denoted by /j TO . Note 
that either scheme can still be used with /j above or below the matching point, it is just that the equations 
(19) are only enforced at /j m . In principle, /j m can be chosen at any value - different choices lead to the same 
overall results, up to higher order corrections. As can be seen in Eqs. 20 and 21, in the generalized MS scheme, 
Aa(/i) and Af a (x,/i) are both of the form a s (/i)\n(/i/m c ) C\ + 0(a 2 s ). Thus, if one chooses /j TO = m c , both 
functions and f a (x,n) in the 3-flavor scheme are equal to their counterpart in the 4-flavor scheme to first 
order in a s at the matching point. Most recent works adopt this choice [9, 41, 18]; we do the same in this paper. 
Although this choice is convenient, it is not required in the general formalism. The ideas behind the matching 
conditions, Eq. 19, are illustrated in Fig. 8 which show two possible matching points, the first being the special 
one [i m = m c . This plot also shows that /j m should not be chosen too far above m c , lest the factor a s ln(/j/m c ) 
in the discontinuity ceases to be perturbative. 

• Choice of transition scale: With F^ X (Q 2 ,x,m c , fi(x,Q,m c )), i = 3,4, we have two sets of calculations, one 
for each scheme. For physics applications, we need to specify which of these to use, say 

j?phys, n 2 r m x_ r F^ x (Q 2 ,x,m c ,n(x,Q,m c )) Q < Ht . , 

*a,\ W ,x,m c )-\ (4) (ZZ) 
F A,\(Q ,x,m c ,fi(x,Q,m c )) Q > [i t 

where we have introduced another scale Ht - the transition point - where one switches from one scheme to the 
other, according to which one is more appropriate, as discussed in the introduction, cf. Fig.l and the previous 
subsections of this appendix. 

Conceptually, the transition point is distinct from the matching point, as should be clear from their defining 
equations, 19 and 22.* The guiding principle for choosing the optimal Ht is that it should be within the region 
where the Ff \ and F^\ calculations are both valid, and that their differences within this region are small. For 
instance, in the idealized situation depicted in Fig. 1, fi t is best chosen to be in the middle of the Q range where 
both uncertainty bands are relatively narrow. 

In practice, one can only estimate the range of uncertainties of the 3-flavor and 4-flavor calculations, say by 
examining the scale-dependence of the respective calculations and then make a judicious choice of fi t - In the 
case of inclusive charm production discussed in Sec. 3, Fig. 5 shows that the transition scale can be chosen at a 
relatively low value, close to m c . For this choice, the composite scheme calculation reduces, in practice, to just 
the 4-flavor calculation. 



oc s (n) or f(x,\x) 

4-flavor 




A ~ a s log^mj + a 2 



Figure 8: Schematic illustration of the matching between 3-flavor and 4-flavor schemes. Two possible matching 
points are shown: the first one, /i m i , is at the the value m c where the discontinuity is of order a 2 ; the second one, 
/i m 2, is at an arbitrary value where the discontinuities are of order a s ln(p/m c ). 

A. 5 What do "LO" and "NLO" mean? 

As already mentioned in the body of this paper, in a multi-scale problem such as heavy quark production, the 
designation of "LO" and "NLO" to a given calculation can be rather misleading in conventional fixed-order calcu- 
lations, due to the presence of large logarithms which vitiates the naive counting of powers of a s . In the composite 
scheme, which has a wider range of applicability than FFN schemes, the meaning of "LO" and "NLO" can be 
better defined, provided the relative magnitudes of the large scales are properly kept in mind. We elaborate a little 
bit. 

• In the 3-flavor scheme, "LO" consists of the 0{a\) HC -y*g — > cc process; whereas "NLO" involves 0(a 2 ) 
processes such as -y*g — > gcc. This formal designation makes physical sense only in the threshold region. 

• In the 4-flavor scheme, the "LO" process is represented by the O(a s ) HE 7*c — > c process; and the "NLO" ones 
consist of the 0(a\) HC process as well as the 0(a\) HE -y*c — > gc process. This formal designation coincides 
with the familiar one in the conventional treatment of DIS structure functions; it makes physical sense only when 
Q 2 > m 2 . 

The apparent mismatch of orders in a s (e.g. order a s being "LO" in the former, but "NLO" in the latter) can be 
understood within the composite scheme which takes into account the order of magnitudes of the other relevant 
quantities in the factorization formula at the appropriate energy ranges. 

Specifically, the 0{a\) and 0(a 2 ) 3-flavor calculations lose their "LO" and "NLO" meaning as Q 2 becomes very 
large compared to m 2 , since each power of a s is accompanied (and neutralized) by a large logarithm \n(Q 2 /m 2 ) 
factor in the hard cross-section. Consequently, all terms become effectively O(a s )\ In fact, the resummation of 
these large logarithms to all orders of a s gives rise precisely to the charm distribution, resulting in the 0(a®) 
charm excitation process 7*c — > c of the 4-flavor scheme. Conversely, in the 4-flavor scheme, although the charm 
quark distribution f%{x, Q) can be considered to be O(a s ) at high Q 2 (where all flavors are on an equal footing), 
as we go down to a lower energy range Q 2 ~ m' 2 , one finds f%{x,Q) ~ a\ ln(Q 2 /m 2 ) ~ 0(a*) compared to 
the dominant gluon and light quark distribution functions (assuming there is no large non-perturbative charm 
component). Therefore, to consistently match the 4-flavor calculation onto the LO (0(a*)) 3-flavor calculation, 
one must include both the "LO" O{a s ) 7*c — > c and the "NLO" 0(a\) j*g — > cc contributions, along with the 
associated subtraction term in the 4-flavor calculation (Cf. Ref. [12]). This implies that our 4-flavor calculation, 
which is NLO at high energy scales, becomes effectively "LO" near the threshold, because both the formally O(a s ) 
and 0(a\) contributions become of the same order of magnitude, 0(a).) - as in the LO 3-flavor scheme. (The 
calculations in the main part of this paper show that, with the natural choice of scale /j, = \JQ 2 + m 2 , the numerical 
predictions of this calculation are actually fairly close to those of the NLO 3-flavor calculation, which is consistent 
with this observation because the size of the NLO correction is within the range of uncertainty of a LO result.) 

This mixing of terms with different apparent powers of a s is physically natural ( cf. Fig. 1 ) and logically 
consistent - it is a necessary feature of switching between different primary schemes, since any finite renormalization 

alwavs pnt.ails a rpsiim.m.n.H.nn. i\ p rp-nro"ani7at.inn , l nf flip nprt.nrhat.inn sprips tn all nrrlprs Tn mnrp rnnrrpt.p t.prms 



the need for mixing terms of different apparent powers of a s arises when: 

(i) the LO diagrams for different subprocesses start at different orders of a s ; 

(ii) the associated parton densities are of different numerical orders of magnitude (such as between g, q, c); 

(iii) the order of magnitude of a parton distribution changes as it evolves with Q (such as for c in the region above 
the threshold); and 

(iv) the hard cross-section contains logarithms of ratios of energy scales which become large." 

The generalized MS formalism, by keeping the physical m c , provides the appropriate scheme to describe the 
underlying physical processes in the different regions encountered in heavy quark production. 
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